一边学习，一边总结，一边分享！
由于微信改版，一直有同学反映。存在长时间接收不到公众号的推文。那么请跟随以下步骤，将小杜的生信筆記设置为星标，不错过每一条推文教程。
欢迎关注《小杜的生信笔记》！！
如何加入社群
小杜的生信笔记，仅有微信社群。
1. 微信群：付费社群。添加小杜好友，加友请知：加友须知！！，加入社群请查看小杜生笔记付费加友入群声明。
2. 小杜个人微信：若你有好的教程或想法，可添加小杜个人微信。值得注意的是，小杜个人微信并不支持免费咨询长时间咨询，但支持小问题2-3个免费咨询。
小杜微信：
知识星球：
2022年教程总汇
2023年教程总汇
11. 靶基因预测与定量
11.1 植物靶基因预测
使用psRobot和targetfinder两款软件
基于targetfinder软件
mamda install targetfinder
## 使用targetfinder进行靶基因预测
targetfinder_threads.pl \
-f all.mature.fa \ # miRNA文件
-d exon.fa \ # mRNA序列文件
-o targetfinder.out \ # 输出结果
-p table # 输出结果格式
基于psRobot
# 软件下载
wget http://omicslab.genetics.ac.cn/psRobot/program/WebServer/psRobot_v1.2.tar.gz
# 解压
tar -zxvf psRobot_v1.2.tar.gz
# 安装软件
./configure
make
make install
## 使用psRobot进行植物靶基因预测
/pub/software/psRobot_v1.2/bin/psRobot_tar \
-s all.mature.fa \ # miRNA文件
-t exon.fa \ # mRNA文件
-o psRobot.out # 输出结果
12. miRDeep2定量、miRDeep-P2 miRNA进行预测
12.1 miRDeep-P2软件进行植物miRNA预测。
## 基于small rna 分类注释信息进行提取
$ awk '$NF=="unknown" || $NF == "intron_sense" \
|| $NF=="intron_antisense" ' \
../S2.reads_mapping/P11.out.read_anno.txt > P11.reads.anno.txt
## 将提取结果转成fa格式
$ awk '{print ">"$1"\n"$2}' P11.reads.anno.txt > P11.reads.fa
ls*.reads.anno.txt | while read line; do \
file_name=$(basename "$line");prefix="${file_name%%.reads.anno.txt}" \
awk '{print ">"$1"\n"$2}' "$line" > "${prefix}.reads.fa";done
运行miRDP2程序进行新miRNA预测
bash /pub/software/miRDP2-v1.1.4/1.1.4/miRDP2-v1.1.4_pipeline.bash \
--genome ../S0.ref_prepare/03.genome/genome.fa \ # 基因组文件
--index ../S0.ref_prepare/03.genome/genome \ # bowtie index
--fasta \ # 输入数据为fasta格式
--input ./B01.mapped.fa \ # 数据输入数据名称
--output ./ # 输出结果目录
输出结果为：
B01.mapped_filter_P_prediction ：新 miRNA 预测结果
B01.mapped_filter_P_prediction.bed ：预测结果 bed 格式文件
将每个样本预测到的 novel miRNA进行合并，根据 第七列 miRNA 是否相同进行去重。
# 这个做法不是很合适，有更好的方法，请更正
cat *reads/*filter_P_prediction > all_reads_filter_P_prediction
## 修改预测的新miRNA的名字
perl ../script/novel_name_mirdp2.pl \
-pbed all_reads_filter_P_prediction \ #miRDP2预测结果，并经过过滤重复后文件
-outpre novel_out # 数据结果前缀
12.2 miRDeep2进行miRNA 、novel miRNA定量
# 合并 已知成熟miRNA、novel miRNA与前体miRNA 的序列
cat $kn_mat $nv_mat > all.mature.fa
cat $kn_hp $nv_hp > all.hairpin.fa
## 将所有样品的srna数据合并
cat ../S1.reads_filter/uniq_data/P1*.fa >all.reads.fa
## 使用mirdeep2软件包中的quantifier.pl程序进行表达定量
/pub/software/mirdeep2-0.1.3/bin/quantifier.pl \
-p all.hairpin.fa \ # 前体序列
-m all.mature.fa \ # 成熟体序列
-r all.reads.fa \ # 输入的srna数据文件，fa格式
-g 0 # bowtie比对错配数
12.3 RdDM进行定量时使用
#!/bin/bash
ShortStack --genomefile /home/yanghj/my_data/worksapce/SmallRNA/Small/data/genome/genome.fa \
--known_miRNAs /home/yanghj/my_data/worksapce/SmallRNA/Small/S0.ref_prepare/01.known_miRNA/can.mature.fa \
--threads 30 --dicermin 20 --dicermax 24 --outdir can_result \
--readfile \
/home/yanghj/my_data/worksapce/SmallRNA/Small/S1.reads_filter/clean_data2/PBC688_inc_1_trimmed.fq \
/home/yanghj/my_data/worksapce/SmallRNA/Small/S1.reads_filter/clean_data2/PBC688_inc_2_trimmed.fq \
/home/yanghj/my_data/worksapce/SmallRNA/Small/S1.reads_filter/clean_data2/PBC688_mock_1_trimmed.fq \
/home/yanghj/my_data/worksapce/SmallRNA/Small/S1.reads_filter/clean_data2/PBC688_mock_2_trimmed.fq \
/home/yanghj/my_data/worksapce/SmallRNA/Small/S1.reads_filter/clean_data2/st_8_inc_1_trimmed.fq \
/home/yanghj/my_data/worksapce/SmallRNA/Small/S1.reads_filter/clean_data2/st_8_inc_2_trimmed.fq \
/home/yanghj/my_data/worksapce/SmallRNA/Small/S1.reads_filter/clean_data2/st_8_mock_1_trimmed.fq \
/home/yanghj/my_data/worksapce/SmallRNA/Small/S1.reads_filter/clean_data2/st_8_mock_2_trimmed.fq
13. 靶基因的预测
提取 20-24nt的碱基，进行靶基因预测。
# 注意将其中的 U 碱基替换成 T 碱基 方便后续进行靶基因的对应，也可以不转化，软件能识别 RNA序列
awk 'NR > 1 && $20 >= 20 && $20 <= 24 {print ">"$2"\n"$11}' Results.txt | sed 's/U/T/g' > mature_20_24.fa
# 查看碱基 数目 只保留20-24的
less -S Results.txt|awk ' NR>1 {print ">"$2"\n"$11}'|seqkit fx2tab -l |less
targetfinder进行靶基因预测
# 将提取到的 mature序列，进行靶基因预测 植物靶基因预测可以直接使用 DNA\RNA 序列做输入。
targetfinder_threads.pl -f mature_20_24.fa -d /home/yanghj/my_data/worksapce/SmallRNA/Small/S0.ref_prepare/03.genome/exon.fa -p table -t 30 -o targetfinder.table
psRobot_tar 进行靶基因预测
# 需要先将 mature.fa 转化为tsv,
seqkit fx2tab mature_20_24.fa > mature_20_24.tab
/home/yanghj/my_data/software/psRobot_v1.2/psRobot_tar -s mature_20_24.tab -t /home/yanghj/my_data/worksapce/SmallRNA/Small/S0.ref_prepare/03.genome/exon.fa -o psRobot.out
# 合并俩个预测的结果，将miRNA转化gene 进行靶基因的整理为target2gene
grep -v 'No' targetfinder.table|awk '{print $1"\t"$2}'| sed "s/\.1//g"| sort | uniq > target2gene.txt
grep '>' psRobot.out |sed 's/>//g;s/\.1//g' |awk '{print $1"\t"$4}'| sort | uniq > psRobot2gene.txt
cat psRobot2gene.txt target2gene.txt | sort -k1,1V |uniq > merge_target2gene.txt
# -V 选项来进行“自然排序”（也称为版本排序
# uniq 命令只能去除连续的重复行,得先 sort排序再去重
13. 进行RPM计算、差异分析
## 导入R进行RPM的计算
#cpm标准化，因为miRNA一个reads就是一个count，所以CPM就是RPM
library(tidyverse)
library(readr)
Counts <- read_delim("E:/ShortStack/bw_20_24/Counts.txt",
delim = "\t", escape_double = FALSE,
trim_ws = TRUE) %>%
select(-1,-3) %>%
column_to_rownames("Name")
library(edgeR)
#cpm标准化，因为miRNA一个reads就是一个count，所以CPM就是RPM
Counts_cpm <- cpm(Counts) %>%
as.data.frame() %>%
rownames_to_column("Name")
# 差异分析
perl /home/yanghj/my_data/software/trinityrnaseq-v2.15.1/Analysis/DifferentialExpression/run_DE_analysis.pl \
--matrix ../31.merge/gene.counts.matrix --method DESeq2 \
--samples_file sample_info.txt \
--contrasts contrast.txt
write.table(Counts_cpm, "E:/ShortStack/bw_20_24/Counts_cpm.txt",sep="\t",
row.names = F,col.names = T,quote=F)
# 同样的替换 t(t(Counts)/colSums(Counts) * 1000000)
###--------------------------
# 参考来源 https://www.notion.so/Small-RNA-ac645653ce71488bb33a4b21d363944f?pvs=4#b7dd0dbb7ea04289b0b933bb76390768
<<<>>>
若我们的教程对你有所帮助，请点赞+收藏+转发，大家的支持是我们更新的动力！！
往期部分文章
1. 最全WGCNA教程（替换数据即可出全部结果与图形）
推荐大家购买最新的教程，若是已经购买以前WGNCA教程的同学，可以在对应教程留言，即可获得最新的教程。（注：此教程也仅基于自己理解，不仅局限于此，难免有不恰当地方，请结合自己需求，进行改动。）
2. 精美图形绘制教程
3. 转录组分析教程
4. 转录组下游分析
小杜的生信筆記 ，主要发表或收录生物信息学教程，以及基于R分析和可视化（包括数据分析，图形绘制等）；分享感兴趣的文献和学习资料!!
