第一個任務：任意物種拿到轉錄本ID 與基因ID的對應關系列表
從ensemble數(shù)據(jù)庫中下載了Human的基因組注釋文件gtf文件
1.下載Human的GTF文件
wget -c ftp://ftp.ensembl.org/pub/release-97/gtf/homo_sapiens/Homo_sapiens.GRCh38.97.gtf.gz
2.獲取轉錄本ID和基因ID的對應關系列表
zcat gencode.v29.annotation.gtf.gz|awk '$3=="transcript" {print $o}'|awk -F '\t' '{print $9}'|cut -d';' -f 1,2 |sed -e 's/gene_id//g' -e's/transcript_id//g' -e 's/"http://g' -e s'/;/\t/g' -e s'/^\s//g' > hg38_geneid_tranid_ly.txt
less hg38_geneid_tranid_ly.txt
1564569665.png
zcat gencode.v29.annotation.gtf.gz |awk '$3 =="transcript" {print $0}'|cut -f 9|cut -d';' -f 1,2 | less -S
perl -p -i -e 's/ge.*"(E.*)".*"(E.*)"/$1\t$2/g' 1.tmp
less 1.tmp
1564569666.png
第一種使用了純shell編程，重復使用了sed命令，
第二種使用了shell+perl的方法卑硫。
可以發(fā)現(xiàn)學習編程的好處就是在于可以使用不同的方式去完成相同的目的豁鲤。
看完曾老師的代碼發(fā)現(xiàn)更加簡潔轧铁，所以我準備課后認真學習一下perl命令行乎赴，爭取更加簡便的寫代碼摔竿！
今天在批量獲取轉錄本ID和基因ID的過程種间坐，遇到了一點兒錯誤
1564626654(1).png
發(fā)現(xiàn)第二列是應該是gtf文件種的版本信息吧灾挨。
因此我就查看了這個物種的gtf注釋信息的第9列到底是什么情況。
1564627021(1).png
感覺還是很有意思的竹宋，在這次實操種劳澄，我發(fā)現(xiàn)了物種的gtf注釋文件的格式并不是一成不變的。在gtf框架是一致的情況下蜈七，還有非常細小的差別秒拔。（比如可以增加了版本信息等）
這樣代碼就需要做出一點兒改變:cut -f 1,3 就可以解決了
zcat Bison_bison_bison.Bison_UMD1.0.97.gtf.gz|awk '$3=="transcript" {print $o}'|awk -F '\t' '{print $9}'|cut -d';' -f 1,3 |sed -e 's/gene_id//g' -e's/transcript_id//g' -e 's/"http://g' -e s'/;/\t/g' -e s'/^\s//g' > Bison_bison_bison.Bison_UMD1.0.97.gtf.gz_geneid_tranid_ly.txt
因此給我們的提示就是，在做轉錄本ID和基因ID的對應關系的時候飒硅，要認真查看物種的gtf文件的格式砂缩，那一點兒細微的差別，可能導致批量操作后的結果狡相，并不是你想要的梯轻！
