mirna
是生物中非常重要的一類非編碼小
rna，其在生物體的調控中具有非常重要的作用，在人中大約三分之一的基因受到
mirna
的調控。對於
mirna
轉錄後調控的分析也越來越多。那麼拿到一組
mirna
測序的資料之後我們進行怎樣的分析呢？
第一，對於所有的測序資料，我們都要進行質量的檢測，這裡我常用的檢測軟體是fastqc，（
fastqcdata1.fq -o data1
），得到的結果是乙個資料夾的壓縮形式，裡面可以得到以下所示的資訊：
同時這些結果都有單獨的格式，用於資料展示與質量評定。在新版本的
fastqc
中有乙個新的功能就是識別
reads
中包含的
adapter
序列，並且
fastqqc
中有乙個
adapter
的庫，可以從裡面找到對應的
adapter
序列，再也不用擔心沒有測序報告沒法去掉
adapter了。
第二，對資料進行過濾，這裡我用的是cutadapt，這個軟體可以去掉reads中的adapter，低質量的reads以及過長過短的reads，還可以對reads中含有n的進行處理。（
cutadapt-a agatcggaagag --quality-base 33 -m 10 -q 20--discard-untrimmed -o trim_data1.fqdata1.fq >cutadpt.info
），這裡--discard-untrimmed是把reads中不含有adapter的reads去掉。
第三，由於分析的是mirna-seq，這裡對clean的reads還要進行一下長度分布的統計，一般就是自己寫指令碼，我用的python。
importsys
mirna_len={}
fori inrange(0,52):
if i.startswith('@') ori.startswith('+'):
continue
length = len(i) - 1
mirna_len[length] += 1
fori inmirna_len:
統計完長度分布之後就是做呈現出來，這裡是用r作圖：
#use allreads
s1=read.table("trim_data1.stat")
s2=read.table("trim_data2.stat")$v2
s3=read.table("trim_data3.stat")$v2
s4=read.table("trim_data4.stat")$v2
s5=read.table("trim_data5.stat")$v2
s6=read.table("trim_data6.stat")$v2
data=cbind(s1, s2, s3, s4, s5, s6)
colnames(data)=c("length","data1","data2","data3","data4","data5","data6")
#normalizeby library size
data$kbt_0=100 * data$data1/sum(data$data1)
data$kbt_1=100 * data$data2/sum(data$data2)
data$kbt_3=100 * data$data3/sum(data$data3)
data$kn6_0=100 * data$data4/sum(data$data4)
data$kn6_1=100 * data$data5/sum(data$data5)
data$kn6_3=100 * data$data6/sum(data$data6)
library(reshape2)
data.melt=melt(data, id="length")
library(ggplot2)
pp+geom_line() +
theme( text =element_text(size=30),
panel.background=element_blank(),
axis.line = element_line(size = 1,colour="black"),
axis.text =element_text(colour="black")) +
labs(title="all readslengthdistribution",x="read length", y="fraction (%)")
得到的示意圖如下，一般在21和24nt的位置有兩個峰：
第四，在得到高質量的clean資料之後就是進行比對，將mirna的資料比對到相應物種的基因組上，這裡我用的是bowtie軟體，（
），我分析的植物mirna-seq的資料，比對率超過了90%。
第五，在得到比對的結果之後就是用htseq進行count計數，把在不同的材料中表達的mirna的reads支援數統計出來。
for i indata1 data2 data3 data4 data5 data6 do
htseq-count-s no -t mirna -i id -o $i.hc.sam $i.ht.sammirna_reference/mirna.gff3 | tee$i.count &
ls$i.count>> count.list
done
第六，然後就是重頭戲，差異表達的mirna，這裡分為有重複的處理和沒有重複的處理兩種，對於
沒有重複的用degseq處理，
有重複的用
deseq處理。
沒有重複的用degseq： r
library("degseq")
#bt_0_1
geneexpmatrix1
geneexpmatrix2
write.table(geneexpmatrix1[30:31,],row.names=false)
write.table(geneexpmatrix2[30:31,],row.names=false)
pdf(file="data1_2.pdf")
layout(matrix(c(1,2,3,4,5,6),3, 2, byrow=true))
par(mar=c(2,2, 2, 2))
degexp(geneexpmatrix1=geneexpmatrix1,genecol1=1, expcol1=2,grouplabel1="data1",
geneexpmatrix2=geneexpmatrix2,genecol2=1, expcol2=2,grouplabel2="data2",
method="mars",outputdir="05demirna/degseq")
dev.off()
有重複的用deseq：r
library("deseq")
data=read.table("ht.genotype_data.txt",header=true,row.names=1)
pd=data.frame(row.names=colnames(data),condition=c("data3","data3","data4","data4"),libtype=c("single-end","single-end","single-end","single-end"))
ps=pd$libtype=="single-end"
ct=data[,ps]
condition=pd$condition[ps]
cds=newcountdataset(ct,condition)
cds=estimatesizefactors(cds)
sizefactors(cds)
cds=estimatedispersions(cds)
res=nbinomtest(cds,"data3","data4")
write.table(res,file="data3_data4.xls")
quit()
第七，
在得到的結果中兩種軟體**到的共同的部分是結果比較可信的部分。
第八，當然對於mirna還有很過方面，對靶基因的功能分析，對mirna二級結構的分析，對樣品中新mirna的分析等等。
使用mirDeep2進行miRNA seq資料分析
git clone mirdeep2.0.1.2 cd mirdeep2.0.1.2 使用install.pl指令碼進行安裝 perl install.pl會有如下的提示資訊 提示資訊 可以按照他的要求，直接使用source bashrc載入環境變數，然後再次執行perl install.pl就會幫...
kernel power off流程分析
凡是 linux 核心上層關機時，底層均會調到 kernel power off 電腦可以使用按鍵 ctr alt del 鍵進入關機，下面我們看看 流程 syscall define4 kernel power off pm power off prepare machine power off ...
UBoot流程分析
uboot程式分析 程式入口分析 第一階段bl1程式分析 第二階段bl2程式分析 解壓uboot原始碼，開啟頂層makefile，每個uboot所支援的開發板在makefile中都會有乙個配置選項，在e uboot board samsung smdk2440，有乙個uboot.lds鏈結器指令碼檔...