将fasta格式参考基因组的每一条染色体都划分成一个一个等长的区

昨天晚上，有个师弟给我发消息说，说他工作遇到一问题，尝试了一段时间也没有解决。问题其实很简单，就是想把某个参考基因组的每一条染色体都划分成等长的区间，每个区间2kb，这样的程序怎么写。

其实就是一个小脚本的事情，于是我就随手给他写了一个脚本如下。

#!/usr/bin/env python
# coding=utf-8
'''
@Author: xuzhongtian
@LastEditors: Xu Zhongtian
@email: xuzhongtian11@163.com
@github: https://github.com/BiocompZTXu
@Date: 2019-11-16 12:37:33
@LastEditTime: 2019-11-16 13:19:21
@motto: Still water run deep
@Description: scatter genome sequence to bins with assigned length
@FilePath: /Python/home/xuzhongtian/huyiqing/Bining_genome.py
'''


import argparse

def get_parser():
    '''
    Using argpaser to parse the input parameters
    '''
    parser = argparse.ArgumentParser(description='Scattering genome sequence to bins with assigned length')
    parser.add_argument('-i', '--input', type=str, nargs=1,help='input the reference genome file with fasta format')
    parser.add_argument('-b', '--binsize', type=int, nargs=1, help='bin size')
    return parser


def binGenome(fasta, binsize):
    '''Bining the reference to segments head to end with equal length'''
    sequence = {}
    with open(fasta) as fa:
        for line in fa:
            if line.startswith(">"):
                ac = line.strip().split(">")[1]
                seq = ""
            else:
                seq += line.strip()
            sequence[ac] = seq

    for ac, seq in sequence.items():
        seqlen = len(seq)  
        for bin in range(0, seqlen, binsize):
            start, end = bin, bin + binsize
            if bin+binsize > seqlen:
                end = seqlen
            print(">%s_%d-%d"%(ac,start,end))
            print(seq[start:end])

最后脚本的运行方式如下所示（假设我们需要将这个测试文件切割成一段一段长度为20的bin区间）。

python Bining_genome.py  -i test_genome.fa -b 20 > result.fasta

写脚本的时候，我一般会取一个很小很小的文件做测试，当确定运行正确之后，再用真实的数据来跑。

我们可以看看测试的数据，test_genome.fasta长啥样。

>Chr1
GTTTGGTTTGTGCGTGATGTTAAGATCGGAAGAGCACACGTCTGGAGCACACGTCTG
>Chr2
NCCTCTGCAAACGGGTCTGATAGTATTTCAGATCGGAAGAGCACACGTCT

当运行完这个脚本之后，结果如下：

>Chr1_0-20
GTTTGGTTTGTGCGTGATGT
>Chr1_20-40
TAAGATCGGAAGAGCACACG
>Chr1_40-57
TCTGGAGCACACGTCTG
>Chr2_0-20
NCCTCTGCAAACGGGTCTGA
>Chr2_20-40
TAGTATTTCAGATCGGAAGA
>Chr2_40-50
GCACACGTCT

目的达成。

Done.