matplotlib:根据gtf文件画出基因全部的转录本结构

import numpy as np
from matplotlib import pyplot as plt

面向对象的界面:显式创建图形和轴,并在它们上调用方法
x = np.linspace(0, 2, 100)
fig, ax = plt.subplots()
ax.set_xlabel('x label')  # x轴
ax.set_ylabel('y label')  # y轴
ax.set_title("Simple Plot")  # 图名
ax.plot(x, x, label='linear')  # 画x
ax.plot(x, x**2, label='quadratic')  # 画x方
ax.plot(x, x**3, label='cubic')  # 画x立方
ax.legend()  # 画图例

pyplot 界面:依靠pyplot自动创建和管理图形和轴,并使用pyplot函数进行绘图
x = np.linspace(0, 2, 100)
plt.xlabel('x label')
plt.ylabel('y label')
plt.title("Simple Plot")
plt.plot(x, x, label='linear')  # Plot some data on the (implicit) axes.
plt.plot(x, x**2, label='quadratic')  # etc.
plt.plot(x, x**3, label='cubic')
plt.legend()
image.png
from matplotlib import pyplot as plt

fig = plt.figure(1)  # 新建一个名叫 Figure1的画图窗口
fig.patch.set_alpha(1)  # 设置不透明度,默认为1,完全不透明
fig.patch.set_facecolor('w')  # 自定义背景色 "w" 白色

根据gtf文件画出某个基因全部的转录本结构

http://www.biotrainee.com/thread-624-1-1.html 此片段改编自生信技能树论坛生信编程直播第五题:根据GTF画基因的多个转录本结构 发表于 2017-10-4 00:41:45 的 源氏 的回答

# https://ftp.ensembl.org/pub/release-87/gtf/homo_sapiens/Homo_sapiens.GRCh38.87.chr.gtf.gz
def find_target_data(gtf, gene):
    from collections import defaultdict
    gene_transInfo = defaultdict(list)
    with open(gtf, 'r') as f:
        """1    ensembl_havana  gene    65419   71585   .       +       .       gene_id "ENSG00000186092"; 
                gene_version "6"; gene_name "OR4F5"; gene_source "ensembl_havana"; gene_biotype "protein_coding";"""
        for line in f:
            if f'gene_name "{gene.upper()}"' in line:
                line_spt = line.strip().split('\t')
                try:
                    chr, db, record, start, end, score, strand, phase, info = line_spt
                    gene_transInfo['start'].append(start)
                    gene_transInfo['end'].append(end)
                    gene_transInfo['record'].append(record)
                except:
                    pass
    if not gene_transInfo:
        print('\n\n There is some wrong with your gene name!\n')
        raise NameError('your gene is not exit')

    return gene_transInfo


def draw_gene_structure(gene, gene_transInfo,  png_path):
    from matplotlib import pyplot as plt
    from matplotlib.ticker import MultipleLocator
    color_dic = {'exon': '#00896C', 'CDS': '#78552B', 'start_codon': '#CB1B45', 'stop_codon': 'black',
                 'five_prime_utr': '#F19483', 'three_prime_utr': '#2EA9DF'}
    linewith_dic = {'exon': 8.0, 'CDS': 6.0, 'start_codon': 8.0, 'stop_codon': 8.0, 'five_prime_utr': 4.0, 'three_prime_utr': 4.0}

    gene_start = min(map(int, gene_transInfo['start'])) - 500   # 两边边扩大些
    gene_end = max(map(int, gene_transInfo['end'])) + 500

    fig = plt.figure(figsize=(12, 6))  # 建个画板,画板大小
    ax = fig.add_subplot()  # 画板上建个画布,这里可以建多个画布
    ax.set_xlim(int(gene_start) - 500, int(gene_end) + 500)  # 画x轴
    ax.ticklabel_format(useOffset=False, style='plain')  # x,y轴禁用科学记数法

    t = 0
    record = gene_transInfo['record']
    for i in range(len(record)):
        if record[i] == 'transcript':
            t += 1
            ax.plot([int(gene_transInfo['start'][i]), int(gene_transInfo['end'][i])], [t, t], color="black")
        elif record[i] == 'gene':
            pass
        else:
            ax.plot([int(gene_transInfo['start'][i]), int(gene_transInfo['end'][i])], [t, t], color=color_dic[record[i]], linewidth=linewith_dic[record[i]])

    ax.set_title(f"the transcripts of {gene} gene")
    ymajorLocator = MultipleLocator(1)  # y轴主间距设为1,事实上按数量怎么顺眼怎么来
    ax.yaxis.set_major_locator(ymajorLocator)  #
    ax.set_ylim(0, t + 3)  # +3 单纯是让y轴长点

    # 添加图例
    import matplotlib.patches as mpatches
    legend_lt = []
    for region in color_dic:
        legend = mpatches.Patch(color=color_dic[region], label=region, linewidth=linewith_dic[region])
        legend_lt.append(legend)
    ax.legend(handles=legend_lt)
    #  这个图例画的有点丑 可以再调整调整

    fig.savefig(png_path, dpi=150)
    plt.show()


gtf_path = r'D:\database\hg38\Homo_sapiens.GRCh38.100.chr.gtf'
gene = 'ATP7B'
gene_transInfo = find_target_data(gtf_path, gene)
png_path = f'D:\result\{gene}.png'
draw_gene_structure(gene, gene_transInfo, png_path)

image.png

调整坐标轴刻度间隔 (本片段为CSDN博主「南辰以北」的原创文章)

from matplotlib.ticker import MultipleLocator, FormatStrFormatter

xmajorLocator   = MultipleLocator(1)   x轴主刻度
ax.xaxis.set_major_locator(xmajorLocator)

ymajorLocator   = MultipleLocator(1)   y轴主刻度
ax.yaxis.set_major_locator(ymajorLocator)

xminorLocator   = MultipleLocator(0.25)   x轴副刻度
ax.xaxis.set_minor_locator(xminorLocator)

yminorLocator   = MultipleLocator(0.25)   y轴副刻度
ax.yaxis.set_minor_locator(yminorLocator)
————————————————
版权声明:本文为CSDN博主「南辰以北」的原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接及本声明。
原文链接:https://blog.csdn.net/weixin_34498545/article/details/112631706
最后编辑于
?著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 213,992评论 6 493
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 91,212评论 3 388
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 159,535评论 0 349
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 57,197评论 1 287
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 66,310评论 6 386
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 50,383评论 1 292
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 39,409评论 3 412
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 38,191评论 0 269
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 44,621评论 1 306
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 36,910评论 2 328
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 39,084评论 1 342
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 34,763评论 4 337
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 40,403评论 3 322
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 31,083评论 0 21
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 32,318评论 1 267
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 46,946评论 2 365
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 43,967评论 2 351

推荐阅读更多精彩内容