要處理四個矩陣文件遮婶，分別為GSE109555 GSE136447 GSE36552 GSE66507
從ncbi網(wǎng)站上下載原始.sar文件政供，以及轉(zhuǎn)化壓縮成為.fastq.gz文件
在后臺處理利用nohup&或者screen
cat GSE109555.txt | while read i;
do
/data/biosoft/software/sratoolkit.2.9.6-1-ubuntu64/bin/prefetch -X 9999999999999999999 -O `pwd` $i && echo "**${i}.sra done**";
fastq-dump --split-3? ${i}.sra;
gzip?${i}.fastq
done
————————————生成count矩陣渴逻，要求比對率為60%以上————————————
GSE66507，是雙端測序文件精钮，長度為50bp
比對過程中出現(xiàn)問題
多是由于數(shù)據(jù)沒有下載完整造成的袍啡，此時需要重新下載
先是利用組里的index進行比對；結(jié)果還可以钾怔？大致有2/3比對率是大于百分之50的
GSE36552碱呼，是單端測序文件，長度為100bp
準確率較好宗侦，基本可以穩(wěn)定在55%以上
但是對于
GSE109555愚臀，是雙端測序文件，長度為150bp
以及
GSE136447凝垛，是雙端測序文件懊悯，長度為150bp
比對率較差蜓谋，都只有15%左右
更換比對文件Homo_sapiens.GRCh38.99.gtf;hg38_ercc.fa
STAR建立index文件
STAR --runMode genomeGenerate --genomeDir /data/shift/other/index_STAR --runThreadN 20 --genomeFastaFiles hg38_ercc.fa --sjdbGTFfile Homo_sapiens.GRCh38.99.2.gtf
出現(xiàn)問題：.fa和gtf文件里的chr1 和1之間有差異，求助解決
Solution: check the formatting of the GTF file. Most likely cause is the difference in chromosome naming between GTF and FASTA fill
確認文件中的開頭
grep -o -E "^\w+([-+.]\w+)*" Homo_sapiens.GRCh38.99.gtf|sort|uniq
將.gtf文件中加chr
繼續(xù)建立index文件，而后進行比對：本次比對對數(shù)據(jù)進行了清洗
#!/bin/bash
cat /data/shift/other/GSE36552/GSE36552.txt|while read i;
do
fastp -i /data/shift/other/GSE36552/GSE36552_data/${i}.fastq.gz? -o /data/shift/other/GSE36552/GSE36552_data_clean/${i}_clean.fastq.gz;
STAR? --runThreadN? 12? --genomeDir /data/shift/other/index_STAR? --readFilesCommand? zcat --readFilesIn? /data/shift/other/GSE36552/GSE36552_data_clean/${i}_clean.fastq.gz? --sjdbGTFfile /data/shift/other/index_STAR/Homo_sapiens.GRCh38.99.2.gtf --outFileNamePrefix /data/shift/other/GSE36552/GSE36552_count.38.99/${i}_ --outSAMtype BAM SortedByCoordinate;
featureCounts -p -t exon -g gene_id -a /data/shift/other/index_STAR/Homo_sapiens.GRCh38.99.2.gtf -o /data/shift/other/GSE36552/GSE36552_count.38.99/${i}_counts.txt /data/shift/other/GSE36552/GSE36552_count.38.99/${i}_Aligned.sortedByCoord.out.bam;
done
但是比對結(jié)果并不盡如人意捧毛，比對成功率略低于組里文件：思考建立STAR index時是否應(yīng)該根據(jù)不同的序列長度增加參數(shù) --sjdbOverhang?reads長度的最大值減1观堂，默認是100
提出解決辦法：隨機選取一對文件fastq從中取出50bp(兩個文件對稱選取)的長度，與人進行比對
先對GSE136447嘗試
利用軟件seqkit
從Download - SeqKit - Ultrafast FASTA/Q kit (shenwei.me)下載軟件
報錯：seqkit無法處理.gz文件
報錯：更換剪切長度50咐吼，100都報錯
截取后文件分析
比對源文件
原文件就有點問題（心酸）
嘗試解決途徑吹缔，僅保留后50位數(shù)據(jù)，報錯不行
更改保留前100數(shù)據(jù)嘗試,比對成功
突然發(fā)現(xiàn)小時錯了/(ㄒoㄒ)/~~看比對率的文件應(yīng)該是锯茄，*_Log.final.out
GSE66507的比對率為80%以上
GSE36552的比對率為87%以上
現(xiàn)有問題GSE136447 featureCounts的assignment低
按照文章方式進行比對
利用hisat2
我知道為啥我的GSE109555一直都不好使了
原因一：有好多數(shù)據(jù)都沒下完整（可能是由于網(wǎng)絡(luò)原因）厢塘，這部分需要重新下載
原因二：文章中說不能直接用star，需要將每一個文件拆分后肌幽，將子文件做star
1.利用 zcat SRR8569907_2.fastq.gz | wc -l 與原文件比對（一個一個比啊晚碾，愿天堂沒有數(shù)據(jù)不完整）
2.下載 GSE109555_TrioSeq_DataInfo.txt.gz 作為輸入執(zhí)行
perl s01.Barcode_UMI_QC_per1w_V2.pl test_1.fq.gz test_2.fq.gz 不同的文件對應(yīng)不同應(yīng)當根據(jù)SRR對應(yīng)的GME查找 ./
打開文件夾以后里面對應(yīng)的是每個壓縮文件
而后就可以執(zhí)行star了！
85%還挺好牍颈，可以說是感動中國了/(ㄒoㄒ)/~~接下來就是等待~
