最近在自学linux，看到生信技能树分享的linux习题（http://www.bio-info-trainee.com/2900.html）就稍微动手做了一下，感觉要学习的地方还有很多。在此记录一下自己的答案。
一、在任意文件夹下面创建形如 1/2/3/4/5/6/7/8/9 格式的文件夹系列。
$ mkdir -p 1/2/3/4/5/6/7/8/9
二、在创建好的文件夹下面，比如我的是 /Users/jimmy/tmp/1/2/3/4/5/6/7/8/9 ，里面创建文本文件 me.txt
$ touch 1/2/3/4/5/6/7/8/9/me.txt
三、在文本文件 me.txt 里面输入内容:
Go to: http://www.biotrainee.com/
I love bioinfomatics.
And you ?
vim 1/2/3/4/5/6/7/8/9/me.txt
#在vim中按i进入INSERT模式，输入内容，`:wq`保存并退出
四、删除上面创建的文件夹 1/2/3/4/5/6/7/8/9 及文本文件 me.txt
$ rm -rf 1
$ ls
五、在任意文件夹下面创建 folder1~5这5个文件夹，然后每个文件夹下面继续创建 folder1~5这5个文件夹，效果如下：
$ mkdir -p folder{1..5}/folder{1..5}
$ tree
六、在第五题创建的每一个文件夹下面都 创建第二题文本文件 me.txt ，内容也要一样。(这个题目难度超纲，建议一个月后再回过头来做)
# shell scipt
#!/bin/bash
#
for i in {1..5};do
for j in {1..5};do
touch folder${i}/folder${j}/me.txt
done
done
七，再次删除掉前面几个步骤建立的文件夹及文件
$ rm -rf folder{1..5}
八、下载 http://www.biotrainee.com/jmzeng/igv/test.bed 文件，后在里面选择含有 H3K4me3 的那一行是第几行，该文件总共有几行。
$ wget http://www.biotrainee.com/jmzeng/igv/test.bed
$ grep -n H3K4me3 test.bed #`-n`显示行号
$ wc test.bed
10 88 3099 test.bed
九、下载 http://www.biotrainee.com/jmzeng/rmDuplicate.zip 文件，并且解压，查看里面的文件夹结构
$ wget http://www.biotrainee.com/jmzeng/rmDuplicate.zip
$ unzip rmDuplicate.zip
$ tree rmDuplicate
十、打开第九题解压的文件，进入 rmDuplicate/samtools/single 文件夹里面，查看后缀为 .sam 的文件，搞清楚 生物信息学里面的SAM/BAM 定义是什么。
$ cd rmDuplicate/samtools/single/
$ less tmp.sam
SAM文件
SAM(Sequence Alignment/Map)格式是一种通用的比对格式，用来存储reads到参考序列的比对信息。
SAM是一种序列比对格式标准，由sanger制定，是以TAB为分割符的文本格式。主要应用于测序序列mapping到基因组上的结果表示，当然也可以表示任意的多重比对结果。
SAM分为两部分，注释信息（header section）和比对结果部分（alignment section）。行：除注释外，每一行是一个read。
BAM文件
BAM是SAM的二进制格式，因此两者格式相同，只是BAM文件占用储存空间更小，运算更快
作者：井底蛙蛙呱呱呱
链接：https://www.jianshu.com/p/9c99e09630da
来源：简书
著作权归作者所有。商业转载请联系作者获得授权，非商业转载请注明出处。
十一、安装 samtools 软件
$ conda install -c bioconda samtools
十二、打开 后缀为BAM 的文件，找到产生该文件的命令。
在@PG这行找到相应产生命令
十三题、根据上面的命令，找到我使用的参
考基因组 /home/jianmingzeng/reference/index/bowtie/hg38 具体有多少条染色体。
hg38有许多组装版本的chromosomes，都不是主要的，所以在结果中将其排除
$ samtools view -h tmp.sorted.bam |grep @SQ | cut -f 2 | sort |uniq | grep -v '_'
SN:chr1
SN:chr10
SN:chr11
SN:chr12
SN:chr13
SN:chr14
SN:chr15
SN:chr16
SN:chr17
SN:chr18
SN:chr19
SN:chr2
SN:chr20
SN:chr21
SN:chr22
SN:chr3
SN:chr4
SN:chr5
SN:chr6
SN:chr7
SN:chr8
SN:chr9
SN:chrM
SN:chrX
SN:chrY
十四题、上面的后缀为BAM 的文件的第二列，只有 0 和 16 两个数字，用 cut/sort/uniq等命令统计它们的个数。
先看看.bam文件中出现“0/16”的行的特征
# ~/practice/biotrainee/data/rmDuplicate/samtools/single
$ samtools view -h tmp.sorted.bam |less
都以SRR开头(但其实只要在samtools view命令后去掉参数-h就只会显示reads信息了)
# ~/practice/biotrainee/data/rmDuplicate/samtools/single
$ samtools view -h tmp.sorted.bam |grep SRR | cut -f 2 | uniq -c
29 0
24 16
十五题、重新打开 rmDuplicate/samtools/paired 文件夹下面的后缀为BAM 的文件，再次查看第二列，并且统计
sam/bam文件中第二列为FLAG，其意义如下图所示
# ~/practice/biotrainee/data/rmDuplicate/samtools/paired
$ 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
实际.bam文件中出现的FLAG并不完全是如上图的整数，而是组合而来的，如147=128+16+2+1
十六题、下载 http://www.biotrainee.com/jmzeng/sickle/sickle-results.zip 文件，并且解压，查看里面的文件夹结构， 这个文件有2.3M，注意留心下载时间及下载速度。
$ wget http://www.biotrainee.com/jmzeng/sickle/sickle-results.zip
$ tree sickle-results
十七题、解压 sickle-results/single_tmp_fastqc.zip 文件，并且进入解压后的文件夹，找到 fastqc_data.txt 文件，并且搜索该文本文件以 >>开头的有多少行？
# ~/practice/biotrainee/data
$ unzip sickle-results/single_tmp_fastqc.zip
# 解压后文件在当前工作文件夹，而非解压包所在文件夹
$ cd sickle-results/
$ cat fastqc_data.txt |grep "^>>"|wc -l
24
十八题、下载 http://www.biotrainee.com/jmzeng/tmp/hg38.tss 文件，去NCBI找到TP53/BRCA1等自己感兴趣的基因对应的 refseq数据库 ID，然后找到它们的hg38.tss 文件的哪一行。
$ wget http://www.biotrainee.com/jmzeng/tmp/hg38.tss
$ cat hg38.tss | grep "NM_002046"
NM_002046 chr12 6532405 6536405 0
GAPDH transcript variant1
十九题、解析hg38.tss 文件，统计每条染色体的基因个数。
$ cat hg38.tss | cut -f 2 | sort | uniq -c| grep -v "_"
6050 chr1
2824 chr10
3449 chr11
2931 chr12
1122 chr13
1883 chr14
2168 chr15
2507 chr16
3309 chr17
873 chr18
3817 chr19
4042 chr2
1676 chr20
868 chr21
1274 chr22
3277 chr3
2250 chr4
2684 chr5
3029 chr6
2720 chr7
2069 chr8
2301 chr9
2 chrM
2553 chrX
414 chrY
这里的思路是统计每条染色体出现的次数来初步估算转录本的个数，如果想统计基因个数还需要知道NM,NR与基因对应关系（因为一个基因可能有多个转录本）。
二十题、解析hg38.tss 文件，统计NM和NR开头的数量（原为“熟练”），了解NM和NR开头的含义。
$ cat hg38.tss | cut -f 1 | cut -d "_" -f 1 | sort | uniq -c
51064 NM
15954 NR
NM: Refseq认证的确凿的mRNA转录本
NR：Refseq认证的RNA转录本，为non-coding RNA
网友评论