Linux系统RNA-seq分析（stringtie转录本组装)

Luminosite Eternelle

已于 2024-03-14 21:12:01 修改

阅读量971

点赞数 23

文章标签： linux 运维

于 2024-03-14 01:00:19 首次发布

本文链接：https://blog.csdn.net/weixin_43927366/article/details/136692728

版权

前言

Stringtie对于转录本组装，如果仅仅关注已知基因的表达情况，不研究新的转录本则可以跳过这步，直接进行定量。

一、转录组组装

#新建trans_ass文件夹，将排序的BAM文件，基因组gtf文件转移到这个文件夹中
#trans_ass文件夹中新建output子文件夹
mkdir -p output  
stringtie -G tair10.gtf -o 006output/SRR3418006-stringtie.gtf SRR3418006-filtered-STAR.sorted.bam

这里，-G tair10.gtf 告诉StringTie使用tair10.gtf作为参考基因组注释文件。-o 005output/SRR3418006.gtf 指定输出文件的名称和位置，这里我们将其放在006output子文件夹中，并命名为SRR3418006.gtf。最后，SRR3418006-filtered-STAR.sorted.bam 是你的输入BAM文件。

二、转录本整合

同上步，将剩余所有其他样本的sorted.bam文件相同的操作得到样本的gtf文件。

echo -e "/home/lumino/TEST1/trans_ass/output/SRR3418005-stringtie.gtf\n/home/lumino/TEST1/trans_ass/output/SRR3418006-stringtie.gtf" > gtflist.txt

echo: 这是一个 shell 命令，用于在终端上显示一个字符串或变量的值。在这里，它被用来生成两个文件路径的字符串。
-e: 这是 echo 命令的一个选项，允许对转义字符进行解释。在这种情况下，虽然命令中没有使用任何转义字符（如 \n 用于换行），但 -e 选项通常用于确保 echo 命令能够正确解释所有可能的转义字符。
"/home/lumino/TEST1/trans_ass/output/SRR3418005-stringtie.gtf\n/home/lumino/TEST1/trans_ass/output/SRR3418006-stringtie.gtf": 这是 echo 命令要显示的字符串。字符串中包含两个文件路径，它们之间用 \n 分隔。\n 是一个转义字符，代表换行符，所以这两个文件路径会在输出的文本中分别位于不同的行。
>: 这是一个重定向操作符。它用于将 echo 命令的输出（即上述字符串）写入一个文件，而不是显示在终端上。如果文件不存在，它会被创建；如果文件已经存在，它的内容会被覆盖。
gtflist.txt: 这是重定向操作符的目标文件名。echo 命令输出的字符串（即两个文件路径）将被写入这个文件。

综上所述，这个命令会创建（或覆盖）一个名为 gtflist.txt 的文件，并在其中写入两行文本，分别是两个 .gtf 文件的路径。这些路径分别指向 SRR3418005-stringtie.gtf 和 SRR3418006-stringtie.gtf 文件。

#确保产生txt，且路径正确
#将所有样本得到的转录本注释信息整合起来
stringtie --merge -G "../tair10.gtf" -o merged.gtf $(cat ./gtflist.txt)

--merge: 这个选项告诉 stringtie 你想要合并多个 .gtf 文件的结果。通常，这些 .gtf 文件是单个样本的转录本组装结果。
-G "../tair10.gtf": 这个选项用于指定一个参考基因组注释文件（在这里是 tair10.gtf）。这个文件包含了基因和转录本的位置和结构信息。在合并过程中，stringtie 会使用这个参考注释来改进和标准化输出的转录本集合。
-o merged.gtf: 这个选项用于指定输出文件的名称和路径。在这里，合并后的结果将被写入 merged.gtf 文件中。
$(cat ./gtflist.txt): 这部分命令是 shell 脚本的一个特性，用于读取 gtflist.txt 文件中的每一行（即每个 .gtf 文件的路径），并将它们作为参数传递给 stringtie 命令。cat ./gtflist.txt 会输出文件的内容，而 $(...) 会捕获这些输出，并将其作为命令行参数。这样，你就无需手动输入每个 .gtf 文件的路径了。

综合起来，这个命令的作用是：使用 stringtie 的 --merge 功能，结合参考基因组注释 tair10.gtf，将 gtflist.txt 文件中列出的所有 .gtf 文件合并成一个统一的 merged.gtf 文件。

三、转录本注释文件的比较

gffcompare -r ../tair10.gtf -o ./gffcp_merged_output merged.gtf

gffcompare：这是命令本身，它调用 gffcompare 工具来执行比较操作。
-r ../tair10.gtf：-r 选项后面跟着的是参考GTF文件的路径。在这个例子中，参考文件是 ../tair10.gtf，它位于当前目录的上一级目录中。这个参考文件通常包含了已知的、经过验证的基因和转录本的注释信息。
-o ./gffcp_merged_output：-o 选项后面是输出目录的路径。在这个例子中，输出目录是 ./gffcp_merged_output，它位于当前目录下。如果 gffcp_merged_output 这个目录不存在，gffcompare 会尝试创建它。所有输出文件（如 .annotated.gtf、.stats、.loci、.tracking 等）都将保存在这个目录中。
merged.gtf：这是要进行比较的预测的GTF文件。这个文件可能包含了由某个转录组组装工具（如Cufflinks、StringTie等）生成的预测转录本。

执行这个命令后，gffcompare 会读取 merged.gtf 文件中的每个预测转录本，并将其与 ../tair10.gtf 文件中的参考转录本进行比较。然后，它会生成一系列输出文件，这些文件可以帮助用户理解预测转录本与参考转录本之间的匹配情况，包括哪些预测转录本是已知的（与参考注释匹配），哪些是新的（不与参考注释匹配），以及它们之间的其他关系。

查看上述 gffcompare 输出的四个文件，您可以使用多种可视化软件，这些软件通常支持 GTF、GFF 等格式，或者可以处理文本数据以生成可视化。以下是一些建议的可视化软件：

1. IGV (Integrative Genomics Viewer)

IGV 是一个强大的基因组浏览器，它支持多种格式的基因组数据，包括 GTF 和 GFF。您可以加载 gffcompare_result.annotated.gtf 文件到 IGV 中，查看预测转录本与参考注释在基因组上的位置关系。

2. UCSC Genome Browser

UCSC Genome Browser 是另一个广泛使用的基因组浏览器，它提供了丰富的基因组数据和可视化工具。虽然它主要设计用于浏览人类基因组，但它也支持其他物种的基因组。您可以将 GTF 文件上传到 UCSC Genome Browser 中进行查看。

3. JBrowse

JBrowse 是一个快速、灵活的基因组浏览器，可以在本地或服务器上运行。它支持 GTF 和其他格式，并且可以通过自定义配置来显示丰富的基因组注释。

4. R 语言和相关的生物信息学包

如果您熟悉 R 语言，您可以使用 GenomicFeatures、Gviz、ChIPseeker 等包来读取和可视化 GTF 文件。这些包提供了丰富的绘图功能，可以定制复杂的基因组可视化。

5. 自定义脚本和工具

如果您需要更高级的可视化或分析功能，您可以编写自定义脚本（如 Python、Perl 等）来处理 gffcompare 的输出文件，并使用像 matplotlib（Python）或 ggplot2（R）这样的绘图库来创建自定义的可视化。

注意事项：

在使用基因组浏览器时，确保您已经加载了相应的参考基因组序列和注释，以便与 gffcompare 的输出进行比较。
对于 gffcompare_result.stats 和 gffcompare_result.tracking 这样的纯文本文件，您可能需要使用文本编辑器或简单的数据分析工具（如 Excel 或 R）来查看和分析它们。
如果您希望在一个集成的环境中进行可视化和分析，可以考虑使用像 Galaxy 这样的生物信息学工作流平台，它提供了许多工具和工作流的集成，可以方便地处理和分析基因组数据。