搜索
查看: 2620|回复: 2

生信编程实战代码注释【第一题:人类基因组的外显子区.....

[复制链接]

10

主题

11

帖子

400

积分

中级会员

Rank: 3Rank: 3

积分
400
发表于 2017-6-30 15:11:09 | 显示全部楼层 |阅读模式
[Python] 纯文本查看 复制代码
'''
首先声明【代码来自生信编程直播第一题:人类基因组的外显子区域到底有多长的python讲师】
我是新手,只是学习了,做了代码注释,如果注释有错,欢迎小伙伴指出,以防祸害别人
本节课目的:提取CCDS.20160908.txt([url]ftp://ftp.ncbi.nlm.nih.gov/pub/CCDS/current_human/CCDS.20160908.txt[/url])文件,[]中外显子的坐标,计算每一个外显子的长度,然后sum
#注释者:默存,20170630
'''
# coding: utf-8
import os
import re
from collections import OrderedDict
#os.chdir("D:")#修改至我们的工作路径

exon_length = 0#初始化
overlapExons = OrderedDict()#定义有序字典{}

with open('CCDS.20160908.txt','r') as f:#打开文档,命名为f,r:读入模式    
    for line in f :#逐行读取 
        if line.startswith('#'):#去掉以#开始的第一行
            continue
        lists=line.rstrip().split('\t')#使用\t将行拆开,只取[]中的数据
        if lists[-2] == '-':#忽略掉lists[-2]中无exon,以-代替的部分
            continue
        lists[-2] =re.sub('\[|\]','',lists[-2])#正则表达式去[,],使用空格分割,\表示转义
        exons =lists[-2] .split(',')#逗号将倒数第二列中的各个exton分开        
        for exon in exons:
            start = int(exon.split('-')[0])#获取一个exon中的起始坐标,int转化为整型
            end = int(exon.split('-')[1])#获取一个exon中的终止坐标,int转化为整型
            coordinate = lists[0]+':'+exon#给每一个exon一个坐标,去多个基因共用重复计算的exon分 
            if coordinate not in overlapExons.keys():#定义一个字典存储exton坐标,如果存储过跳过,否则存进来            
                overlapExons[coordinate]= 1
                #字典overlapExons一个例子,('1: 2100956-2101042', 1)
                exon_length += end-start                
print(exon_length)
#使用python exon_count.py命令执行代码

'''         
这一节补充学习的   
less -S  CCDS.20160908.txt 比较整齐的显示这个txt文件,方便发现数据规律
head -n  CCDS.20160908.txt>tem.txt 取出txt文件中的前n行放到tem.txt中方便发现数据规律
正则表达式(re):[url]https://deerchao.net/tutorials/regex/regex.htm[/url] (正则表达式30分钟入门教程)
路径修改模块(os)
collection模块
python 网页编辑工具jupyter notebook
'''




上一篇:blast2go本地化
下一篇:ExpressionSet类
回复

使用道具 举报

16

主题

61

帖子

515

积分

高级会员

Rank: 4

积分
515
发表于 2017-6-30 15:22:14 | 显示全部楼层
手动点赞!
回复

使用道具 举报

0

主题

1

帖子

53

积分

注册会员

Rank: 2

积分
53
发表于 2017-7-1 21:27:11 | 显示全部楼层
大赞,支持
#os.chdir("D:")#修改至我们的工作路径这个地方应该是手误,多了一个#吧
回复 支持 反对

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

QQ|手机版|小黑屋|生信技能树 ( 粤ICP备15016384号  

GMT+8, 2019-9-16 09:08 , Processed in 0.030526 second(s), 25 queries .

Powered by Discuz! X3.2

© 2001-2013 Comsenz Inc.