Linux基础(二)

紧接linux基础(一)的内容

11、安装 miniconda,建立虚拟环境,在虚拟环境中安装生信必备软件

# miniconda官网下载linux版安装包:Miniconda2-4.7.10-Linux-x86_64.sh*
>./Miniconda2-4.7.10-Linux-x86_64.sh*     #进入安装包所在目录,运行安装包
>conda config --add channels r
>conda config --add channels conda-forge
>conda config --add channels bioconda
>conda create -n rna python=2  #建立python2环境的名叫rna的虚拟conda环境
>souce activate rna    #激活进入rna环境
>conda install -y bwa gatk4 sra-tools fastqc trim-galore star hisat2 bowtie2 subread htseq multiqc samtools salmon     #安装软件
>source deactivate  #退出环境

12、产生bam文件的命令,示例:

>bowtie2-align-s --wrapper basic-0 -p 20 -x /yanke/yanke21/database/index/bowtie2/hg38 -S WT_rep2_Input.sam -U /tmp/41440.unp
>samtools view tmp.sorted.bam | head -n 2     #运行虚拟环境中的软件要先进入conda环境
SRR1042600.42157053     0       chr1    629895  42      51M     *       0       0       ATAACCAATACTACCAATCANTACTCATCATTAATAATCATAATGGCTATA      CCCFFFFFHHHHHJJJJJJJ#4AGHJJIIJJIIIIIJJJJIJIIIIJJIJI     AS:i:-6 XN:i:0  XM:i:2  XO:i:0  XG:i:0  NM:i:2  MD:Z:11C8A30    YT:Z:UU
SRR1042600.42212881     0       chr1    629895  42      51M     *       0       0       ATAACCAATACTACCAATCANTACTCATCATTAATAATCATAATGGCTATA      @@<FDFFBFDHHFJEIIGJI#3AFHGEHEIJIIGIIGGIJIIJIGIIGIIJ     AS:i:-6 XN:i:0  XM:i:2  XO:i:0  XG:i:0  NM:i:2  MD:Z:11C8A30    YT:Z:UU
SRR1042600.12010763     16      chr1    629895  24      51M     *       0       0       ATAACCAATACTTCTAATCAAAACTCATCATTAATAATCATAATGGCTATA      ?4B?1*4DD?11*1*?+22+<3F:3@EC:CC4EA,DEDDDDD?D3B:==+;     AS:i:-10        XN:i:0  XM:i:4  XO:i:0  XG:i:0  NM:i:4  MD:Z:11C0A1C6T29YT:Z:UU
SRR1042600.29629551     16      chr1    629895  40      51M     *       0       0       ATAACCAATACTACCAATCACTACTCATCATTAATAATCATAATGGCTATA      HGF?JJHHFDHHGJJIHDFA+E?JIJJIIHGJJJJJJJHHHHHFFFFFCC@     AS:i:-8 XN:i:0  XM:i:2  XO:i:0  XG:i:0  NM:i:2  MD:Z:11C8A30    YT:Z:UU
SRR1042600.41910745     0       chr1    629896  42      51M     *       0       0       TAACCAATACTACCAATCAANACTCATCATTAATAATCATAATGGCTATAG      CC@FFFFFHHHHGIIHIJJJ#3<CFHCGGIIIJJJJJJJJIGGFHIIJFII     AS:i:-6 XN:i:0  XM:i:2  XO:i:0  XG:i:0  NM:i:2  MD:Z:10C9T30    YT:Z:UU
SRR1042600.14329856     16      chr1    629896  8       18M1I32M        *       0       0       AAACCAAATCCTCCAATCAAATCCTCATCATTAATAATCATAATGGCTATA      #############################@IHHGCE9GHFHHHDDDDD<@@     AS:i:-18        XN:i:0  XM:i:5  XO:i:1  XG:i:1  NM:i:6  MD:Z:0T6T0A2A9A28        YT:Z:UU

13、查看使用的参考基因组 gencode.v29.annotation.gtf.gz 具体有多少条染色体

>zcat gencode.v29.annotation.gtf.gz | grep -v "#" | cut -f 1| sort | uniq -c
250524 chr1
  99843 chr10
 167540 chr11
 159530 chr12
  42767 chr13
  96124 chr14
 100426 chr15
 130878 chr16
 171137 chr17
  48848 chr18
 168227 chr19
 204338 chr2
  61485 chr20
  29652 chr21
  59764 chr22
 172809 chr3
 108328 chr4
 120515 chr5
 124117 chr6
 130671 chr7
  98583 chr8
  97204 chr9
    143 chrM
  89500 chrX
   9064 chrY

14、上面的后缀为BAM 的文件的第二列,只有 0 和 16 两个数字,用 cut/sort/uniq等命令统计它们的个数。

>samtools view tmp.sorted.bam | cut -f 2 | sort | uniq -c 
   29 0
   24 16

15、重新打开 rmDuplicate/samtools/paired 文件夹下面的后缀为BAM 的文件,再次查看第二列,并且统计

>cd ~/rmDuplicate/samtools/paired
>ll
total 52
drwxrwxr-x 2 yanke21 yanke21  4096 Nov 12  2016 ./
drwxrwxr-x 4 yanke21 yanke21  4096 Nov 12  2016 ../
-rw-rw-r-- 1 yanke21 yanke21   354 Nov 12  2016 readme.txt
-rw-rw-r-- 1 yanke21 yanke21  3111 Nov 12  2016 tmp.header
-rw-rw-r-- 1 yanke21 yanke21  4379 Nov 12  2016 tmp.rmdup.bam
-rw-rw-r-- 1 yanke21 yanke21  1886 Nov 12  2016 tmp.rmdup.vcf.gz
-rw-rw-r-- 1 yanke21 yanke21 11649 Nov 12  2016 tmp.sam
-rw-rw-r-- 1 yanke21 yanke21  4543 Nov 12  2016 tmp.sorted.bam
-rw-rw-r-- 1 yanke21 yanke21  1888 Nov 12  2016 tmp.sorted.vcf.gz
>samtools view tmp.sorted.bam | cut -f 2 | sort | uniq -c 
     8 147
     3 163
     1 323
     1 353
     1 371
     1 387
     1 433
     3 83
     2 97
     9 99

16、下载 http://www.biotrainee.com/jmzeng/sickle/sickle-results.zip 文件,并且解压,查看里面的文件夹结构, 这个文件有2.3M,注意留心下载时间及下载速度。

>wget -c http://www.biotrainee.com/jmzeng/sickle/sickle-results.zip 
>unzip sickle-results.zip
>tree sickle-results
sickle-results
├── command.txt
├── single_tmp_fastqc.html
├── single_tmp_fastqc.zip
├── test1_fastqc.html
├── test1_fastqc.zip
├── test2_fastqc.html
├── test2_fastqc.zip
├── trimmed_output_file1_fastqc.html
├── trimmed_output_file1_fastqc.zip
├── trimmed_output_file2_fastqc.html
└── trimmed_output_file2_fastqc.zip

17、解压 sickle-results/single_tmp_fastqc.zip 文件,并且进入解压后的文件夹,找到 fastqc_data.txt 文件,并且搜索该文本文件以 >>开头的有多少行?

>unzip sickle-results/single_tmp_fastqc.zip
>cd single_tmp_fastqc/ && ll
total 296
drwxrwxr-x  4 yanke21 yanke21   4096 Oct  6  2016 ./
drwxr-xr-x 13 yanke21 yanke     4096 Aug 24 11:29 ../
-rw-rw-r--  1 yanke21 yanke21  12048 Oct  6  2016 fastqc_data.txt
-rw-rw-r--  1 yanke21 yanke21   5186 Oct  6  2016 fastqc.fo
-rw-rw-r--  1 yanke21 yanke21 259290 Oct  6  2016 fastqc_report.html
drwxrwxr-x  2 yanke21 yanke21   4096 Oct  6  2016 Icons/
drwxrwxr-x  2 yanke21 yanke21   4096 Oct  6  2016 Images/
-rw-rw-r--  1 yanke21 yanke21    542 Oct  6  2016 summary.txt
>grep -E '^>>' fastqc_data.txt | wc -l
24

18、下载 http://www.biotrainee.com/jmzeng/tmp/hg38.tss文件,去NCBI找到TP53/BRCA1等自己感兴趣的基因对应的 refseq数据库ID,然后找到它们的hg38.tss文件的哪一行。

https://www.ncbi.nlm.nih.gov/gene/7157

>wget -c  http://www.biotrainee.com/jmzeng/tmp/hg38.tss
>less hg38.tss
NM_012279       chr5    177020680       177024680       0
NR_131773       chr5    177020680       177024680       0
NR_131774       chr5    177020680       177024680       0
NM_001308213    chr5    177020680       177024680       0
NM_001308219    chr5    177020680       177024680       0
NM_001318082    chr7    99502644        99506644        0
NM_001291980    chr5    177084872       177088872       0
NM_001283039    chr11   3795819 3799819 0
NM_001283040    chr11   3795819 3799819 0
NM_138964       chr2    68641588        68645588        0
NM_139057       chr15   100340005       100344005       1
>grep 'NM_001143990' hg38.tss
NM_001143990    chr17   7684071 7688071 0
>grep 'NM_001291588' hg38.tss
NM_001291588    chr17   43209827        43213827        0

19、解析hg38.tss 文件,统计每条染色体的基因个数。

cut -f 2 hg38.tss | sort | uniq -c

20、解析hg38.tss 文件,统计NM和NR开头的熟练,了解NM和NR开头的含义

grep  -E '^NM' hg38.tss | wc -l
51064
grep  -E '^NR' hg38.tss | wc -l
15954

NM开头表示转录产物序列,成熟的mRNA序列
NR开头表示非编码的转录子序列,包括RNAs,假基因转子等

©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 203,937评论 6 478
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 85,503评论 2 381
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 150,712评论 0 337
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 54,668评论 1 276
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 63,677评论 5 366
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 48,601评论 1 281
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 37,975评论 3 396
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 36,637评论 0 258
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 40,881评论 1 298
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 35,621评论 2 321
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 37,710评论 1 329
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 33,387评论 4 319
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 38,971评论 3 307
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 29,947评论 0 19
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 31,189评论 1 260
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 44,805评论 2 349
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 42,449评论 2 342

推荐阅读更多精彩内容