這里借助biopython模塊

參考鏈接是 https://biopython.org/wiki/GFF_Parsing

這里BCBio模塊里GFF()函數(shù)解析的內(nèi)容和Bio模塊里SeqIO()函數(shù)解析的內(nèi)容很像

cds和外顯子的關(guān)系

cds 是 coding sequence 的縮寫

具體關(guān)系看下圖 來自鏈接 http://m.itdecent.cn/p/cc5cd7053d6e

image.png

開頭結(jié)尾的外顯子區(qū)可能會(huì)比cds長 ，因?yàn)殚_頭結(jié)尾的外顯子可能包括 UTR，非翻譯區(qū)

處于中間的外顯子和cds等同

首先是根據(jù)gff文件獲取每條染色體的長度

from BCBio import GFF
in_handle = open("tunisia.gff",'r')
for rec in GFF.parse(in_handle):
    for feature in rec.features:
        if feature.type == "region":
            chromo = feature.qualifiers['ID'][0]
            chrid = chromo.split(":")[0]
            chrLen = chromo.split("..")[1]
            print(chrid,chrLen)

這里你的ID可能需要換成其他，這個(gè)得根據(jù)具體gff文件的內(nèi)容定

image.png

獲取gff文件里的基因都有哪些類型

from BCBio import GFF
from collections import Counter

biotype = []

in_handle = open("tunisia.gff",'r')

for rec in GFF.parse(in_handle):
    for feature in rec.features:
        if feature.type == "gene":
            biotype.append(feature.qualifiers['gene_biotype'][0])
            
in_handle.close()
len(biotype)
biotype[0:15]

Counter(biotype)

image.png

統(tǒng)計(jì)每個(gè)蛋白編碼基因有幾個(gè)轉(zhuǎn)錄本

這里需要記住的是每個(gè)feature對(duì)應(yīng)的還有sub_feature這個(gè)是和SeqIO解析genbank文件有差別的地方

gene對(duì)應(yīng)的 sub_features是mRNA,mRNA對(duì)應(yīng)的sub_features是exon或者cds

in_handle = open("pra-2.gff",'r')

for rec in GFF.parse(in_handle):
    for feature in rec.features:
        if feature.type == "gene":
            gene_name = feature.qualifiers["gene"][0]
            i = 0
            for subfeature in feature.sub_features:
                i = i + 1
            print(gene_name, i)

image.png

去除指定基因類型的注釋文件，

比如這個(gè)例子是去除注釋文件中的所有蛋白編碼基因

in_handle = open("tunisia.gff",'r')

fw = open("pra-3.gff",'w')

for rec in GFF.parse(in_handle):
    tmp = rec.features
    i = 0
    index2delete = []
    for feature in rec.features:
        i = i + 1
        if feature.type == "gene" and feature.qualifiers["gene_biotype"][0] == "protein_coding":
            index2delete.append(i-1)
    for counter,index in enumerate(index2delete):
        index = index - counter
        tmp.pop(index)
        
    rec.features = tmp
    GFF.write([rec],fw)
    
fw.close()

這里用到按照索引刪除列表中的元素 ，這個(gè)邏輯暫時(shí)沒有想明白，代碼是

list_given = [1, 2, 3, 4, 5, 6, 7, 8, 9]
index_to_delete = [1, 3, 6]

for counter, index in enumerate(index_to_delete):
    index = index - counter
    list_given.pop(index)
————————————————
版權(quán)聲明：本文為CSDN博主「CAAT9」的原創(chuàng)文章，遵循CC 4.0 BY-SA版權(quán)協(xié)議，轉(zhuǎn)載請(qǐng)附上原文出處鏈接及本聲明。
原文鏈接：https://blog.csdn.net/exm_further/article/details/112251558

好了今天的內(nèi)容暫時(shí)先到這里了

歡迎大家關(guān)注我的公眾號(hào)

小明的數(shù)據(jù)分析筆記本

小明的數(shù)據(jù)分析筆記本 公眾號(hào) 主要分享：1、R語言和python做數(shù)據(jù)分析和數(shù)據(jù)可視化的簡(jiǎn)單小例子；2、園藝植物相關(guān)轉(zhuǎn)錄組學(xué)、基因組學(xué)、群體遺傳學(xué)文獻(xiàn)閱讀筆記；3、生物信息學(xué)入門學(xué)習(xí)資料及自己的學(xué)習(xí)筆記！