用proActiv从RNA-Seq数据中识别活性和替代启动子

约瑟夫·李，德尼兹·德米西奥伊卢，乔纳森Göke

总结

大多数人类基因有多个启动子，控制不同亚型的表达。使用这些替代启动子可以在转录前调节异构体的表达。替代启动子已被发现在广泛的细胞类型和疾病中是重要的。

proActiv是一种从RNA-Seq数据中分析启动子的方法。proActiv使用对齐的读取作为输入，然后为每个注释启动子生成计数和标准化启动子活性估计值。然后，这些估计可以用来确定哪些启动子是活跃的，哪些启动子是不活跃的，以及哪些启动子在不同条件下会改变其活性。

在这里，我们提供了使用proActiv的快速入门指南，以及在两种条件下识别活性启动子和替代启动子的详细工作流程。

如果您在研究中使用proActiv，请注明:

demirciozylu, Deniz等，“泛癌症转录组分析揭示了通过替代启动子的普遍调控。”细胞178.6(2019): 1465-1477。

内容

快速入门:用proActiv量化启动子活性

proActiv通过RNA-Seq数据估计启动子活性。启动子活性定义为每个启动子启动的转录总量。proActiv将BAM文件或连接文件(TopHat2或STAR)和相关基因组的启动子注释对象作为输入。可选参数条件可以提供，描述与每个输入文件对应的条件。在这里，我们演示了proActiv与STAR连接文件(人类基因组GRCh38 Gencode v34)作为输入。由于大小限制，分析仅限于chr1的子集(10,000,000-30,000,000)。

图书馆(高伦雅芙)作为输入的STAR连接文件列表文件< -list.files（执行（“extdata /装饰图案/连接”，包=“高伦雅芙”)，full.names =真正的）描述实验条件的向量条件< -代表（c（“A549”，HepG2的)，每一个=3.）人类基因组GENCODE v34的启动子注释仅限于chr1的子集promoterAnnotation < -promoterAnnotation.gencode.v34.subset结果< -高伦雅芙（文件=文件,promoterAnnotation =promoterAnnotation,条件=条件)

结果是一个summarizeexperiment对象，可以通过如下方式访问:

一个完整的工作流程，以确定替代启动子的使用

在这里，我们提供了一个完整的分步工作流程，用于使用proActiv分析启动子活性，并在不同条件下的样本中识别替代启动子的使用。我们将比较来自2个不同细胞系(A549和HepG2)的样本，以确定替代启动子。

准备输入数据

proActiv使用RNA-Seq数据来量化启动子活性。用户可以选择使用BAM文件或使用制表符分隔的连接文件作为输入，这些连接文件是在使用TopHat2或STAR进行读对齐时生成的。

下面，我们演示跑步高伦雅芙输入STAR连接文件。此数据取自SGNEx项目，并限制在chr1:10,000,000-30,000,000区域。用于比对的参考基因组是Gencode v34 (GRCh38)。这些文件可以在extdata /装饰图案/连接：

文件< -list.files（执行（“extdata /装饰图案/连接”，包=“高伦雅芙”)，full.names =真正的）

启动子注释的制备

为了量化启动子活性，proActiv使用了一组基于基因组注释的启动子。proActiv允许为TxDb对象或GTF文件中的任何基因组创建启动子注释对象preparePromoterAnnotation函数。用户可以选择传递要使用的GTF/GFF或TxDb的文件路径，或者直接使用TxDb对象作为输入。proActiv包括预先计算的人类基因组启动子注释(GENCODE v34)。但是，由于大小限制，注释被限制在chr1:10,000,000-30,000,000区域。用户可以通过下载GTF文件来构建完整的注释GENCODE点击页面并执行以下步骤。

我们演示了使用GTF和TxDb为人类基因组(GENCODE v34)创建受限启动子注释:

##从GTF文件路径gtf。文件< -执行（“extdata /装饰图案/注释/ gencode.v34.annotation.subset.gtf.gz”，包=“高伦雅芙”）promoterAnnotation.gencode.v34.subset<-preparePromoterAnnotation（文件=gtf.file,物种=“Homo_sapiens”）##从TxDb对象txdb。文件< -执行（“extdata /装饰图案/注释/ gencode.v34.annotation.subset.sqlite”，包=“高伦雅芙”）txdb < -loadDb(txdb.file)promoterAnnotation.gencode.v34.subset<-preparePromoterAnnotation（txdb =txdb,物种=“Homo_sapiens”）

的PromoterAnnotation对象有3个槽位:

当为其他物种创建启动子注释时，可以调用genomeStyles要识别物种使用的参数:

的名字(GenomeInfoDb：：genomeStyles（））#>[1]“拟南芥_thaliana”“秀丽隐杆线虫”#>[3]“Canis_familiaris”“Cyanidioschyzon_merolae”"Drosophila_melanogaster" "Homo_sapiens"#>[7]“Mus_musculus”“Oryza_sativa”#>[9]“胡杨/赤毛霉”“褐家鼠”#>[11]“Saccharomyces_cerevisiae”“Zea_mays”

高伦雅芙运行

一旦基因组中的启动子被识别出来，proActiv就会估计每个注释启动子的启动子活性。在这里，我们为GENCODE Release 34加载预先计算的启动子注释。并给出了实验条件高伦雅芙．这个信息允许高伦雅芙总结不同条件下的结果。

promoterAnnotation < -promoterAnnotation.gencode.v34.subset条件< -代表（c（“A549”，HepG2的)，每一个=3.）结果< -高伦雅芙（文件=文件,promoterAnnotation =promoterAnnotation,条件=条件)

结果是一个SummarizedExperiment以原始/归一化启动子计数、绝对/相对启动子活性和基因表达为分析对象:

显示(结果)#>类:总结实验#> dim: 1380# >元数据(0):#> assays(5): promoterCounts normalizedPromoterCounts#> absolutePromoterActivity relativePromoterActivity基因表达#> rownames(1380): 1 2…1379 1380#> rowData names(14): promoterId geneId…A549.class HepG2.class#> colnames(6): sgnex_a549_illumina_replicate1 .run1.子集. sj .out# SGNEx_A549_Illumina_replicate3.run1.subset.SJ.out…# > SGNEx_HepG2_Illumina_replicate4.run1.subset.SJ.out# > SGNEx_HepG2_Illumina_replicate5.run1.subset.SJ.out#> colData names(2): sampleName条件

的rowDataslot按基因存储每个启动子的启动子ID映射和启动子位置(5 '到3 ')。每种条件的平均绝对启动子活性和基因表达也在这里进行了总结。促进者也分为三类。活性< 0.25的启动子被归为非活性启动子，而每个基因中最活跃的启动子被归为主启动子。活跃在较低级别的启动子被归类为次要启动子。

promoterId geneId seqnames 开始 链 internalPromoter promoterPosition txId A549.mean A549.gene.mean HepG2.mean HepG2.gene.mean A549.class HepG2.class
1 ENSG00000000938.13 chr1 27624062 - 真正的 4 ENST0000…… 0.000000 0.00000 0.000000 0.00000 NA NA
2 ENSG00000000938.13 chr1 27626240 - 假 3. ENST0000…… 0.000000 0.00000 0.000000 0.00000 不活跃的 不活跃的
3. ENSG00000000938.13 chr1 27626569 - 假 2 ENST0000…… 0.000000 0.00000 0.000000 0.00000 不活跃的 不活跃的
4 ENSG00000000938.13 chr1 27635185 - 假 1 ENST0000…… 0.000000 0.00000 0.000000 0.00000 不活跃的 不活跃的
5 ENSG00000001460.18 chr1 24379823 - 真正的 7 ENST0000…… 4.721803 20.91493 4.319706 18.81489 NA NA
6 ENSG00000001460.18 chr1 24391679 - 真正的 6 ENST0000…… 4.207194 20.91493 3.999170 18.81489 NA NA

为了更清晰的下游分析，可以去除启动子活性未量化的单外显子转录本。结果可以这样过滤:

通过消除NA启动子计数来去除单外显子转录本/启动子结果< -结果(complete.cases（化验(结果)＄promoterCounts),)

识别替代启动子

为了确定不同启动子的使用情况，proActiv实现了一个线性模型，根据样本条件回归绝对启动子活性和相对启动子活性。启动子绝对活性的变化表明启动子使用的绝对差异。然而，这并不意味着主启动子的替代使用或转换，因为不同条件之间启动子的相对使用可能仍然相同。因此，来自相对启动子活性的信息也必须考虑在内。因此，替代启动子使用的特点是不同条件下启动子绝对和相对活性的显著变化。

候选的替代启动子使用可以通过呼叫来确定getAlternativePromoters．这个函数的参数为结果返回的对象高伦雅芙,referenceCondition，作为与所有其他样品进行比较的参考条件。

alternativePromoters < -getAlternativePromoters（结果=结果,referenceCondition =“A549”）将绝对启动子活性与条件相匹配…拟合相对启动子活性的条件…显示(alternativePromoters)# > upReg美元#> promoterId geneId padjAbs padjRel#> 137 137 ensg00000076864.19 0.004137450 0.03624923#> 138 138 ensg00000076864.19 0.009286492 0.03364840#> 314 314 ensg00000117632.23 0.002433754 0.02968373# ># > downReg美元#> promoterId geneId padjAbs padjRel#> 141 141 ensg00000076864.19 0.007574438 0.0336484

用以下参数检测替代启动子:

分析和可视化替代启动子的使用

在这里，我们提供了上述工作流返回的数据的几种可视化。

替代启动子的使用

为了在不同条件下可视化具有不同启动子使用的基因，功能boxplotPromoters可以使用。该函数输出感兴趣基因的基因表达、绝对启动子和相对启动子活性的箱线图。boxplotPromoters取以下参数:

下面，我们调用boxplotPromotersRAP1GAP:

情节< -boxplotPromoters(因此,“ENSG00000076864.19”）绝对启动子活动的箱线图图书馆(gridExtra)grid.arrange(情节[[1]],情节[[3.]],nrow =1，ncol =2，宽度=c（3.，2）)

的情节对象分别存储绝对启动子活性、相对启动子活性和基因表达在第一、第二和第三插槽的图。由boxplotPromoters反映启动子137和141的替代用法getAlternativePromoters．

启动子类别比例

在这里，我们可视化的分类注释启动子在两个细胞系。在两种细胞系中，类别之间的比例相似，大多数启动子是不活跃的。

图书馆(ggplot2)rdata < -rowData(结果)创建一个总结细胞系和启动子类的长数据框架pdata1 < -data.frame（cellLine =代表（c（“A549”，HepG2的)，每一个=nrow(rdata)),promoterClass =as.factor（c(rdata＄A549.class, rdata＄HepG2.class)))ggplot（na.omit(pdata1))+geom_bar（aes（x =cellLine,填补=promoterClass))+xlab（细胞株的）+ylab（“数”）+实验室（填补=“子类别”）+ggtitle（“推动者的分类”）

按职位划分的主要/次要促销员

主启动子:小启动子比例与启动子位置的对比分析。该分析仅限于具有至少一个活性启动子的多启动子基因。下面，我们生成HepG2细胞系的图。一般情况下，随着启动子位置的增加，主、小启动子比例降低。

因为很多基因都有很多注释启动子，所以我们会折叠启动子##从第5个位置开始，为了简单起见进入一组pdata2 < -as_tibble(rdata)% > %变异（promoterPosition =ifelse(promoterPosition>5，5, promoterPosition))% > %过滤器(HepG2.class%, %c（“主要的”，“小”）)ggplot(pdata2)+geom_bar（aes（x =promoterPosition,填补=as.factor(HepG2.class)),位置=“填满”）+xlab（表达式(启动子～位置～“5”% - > %“3”）)+ylab（“比例”）+实验室（填补=“子类别”）+ggtitle（“HepG2中主要/次要启动子比例”）+scale_y_continuous（打破了=seq（0，1，０．２５)，标签=paste0（seq（0，One hundred.，25)，“%”）)+scale_x_continuous（打破了=seq（1，5)，标签=c（' 1 '，' 2 '，“3”，“4”，“> = 5”）)

主启动子活性与基因表达

主启动子活性和基因表达的比较，通过所有启动子的总和计算。单启动子基因位于对角线上。多启动子基因位于对角线的右侧。下面，我们生成HepG2细胞系的图。该图表明，单一的主启动子通常不能完全解释基因表达，小启动子也有助于基因表达。

获得HepG2的活跃主启动子majorPromoter < -as_tibble(rdata)% > %group_by(geneId)% > %变异（promoterCount =n（））% > %过滤器(HepG2.class= =“主要的”）pdata3 < -data.frame（高伦雅芙=majorPromoter＄HepG2.mean,geneExp =majorPromoter＄HepG2.gene.mean,promoterCount =majorPromoter＄promoterCount)ggplot(pdata3aes（x =geneExp,y =高伦雅芙))+geom_point（aes（颜色=promoterCount),α=0．5）+ggtitle（“HepG2中的主启动子活性与基因表达”）+xlab（“平均基因表达”）+ylab（“平均主启动器活动”）+实验室（颜色=的数量\ n带注释的启动子的）+geom_abline（斜率=1，拦截=0，颜色=“红色”，线型=“冲”）

t-SNE

我们生成一个包含所有激活启动子的t-SNE图。预期，每个细胞系的复制聚集在一起。

图书馆(Rtsne)##删除不活跃的启动子(稀疏行)数据< -化验(结果)＄absolutePromoterActivity% > %过滤器（rowSums(.)>0）数据< -data.frame（t(数据)数据＄示例< -as.factor(条件)set.seed（40）可重复性#tsne.out<-Rtsne（as.matrix（子集(数据、选择=-c(样本))),困惑=1）情节（x =tsne.out＄Y (,1),y =tsne.out＄Y (,2),bg =数据＄样本,asp =1，坳=“黑”，pch =24，cex =4，主要=t-SNE阴谋与推动者\ n至少在一个样本中活跃，xlab =“T-SNE1”，ylab =“T-SNE2”，xlim =c（-300，300)，ylim =c（-300，300）)传说（“topright”，插图=02，title =细胞株的，独特的(条件),pch =c（24，24)，pt.bg =1：长度（独特的(条件)),cex =1．5，电池=“n”）

得到帮助

您可以在Bioconductor支持站点提出疑问和问题:https://support.bioconductor.org．确保你的帖子带有标签高伦雅芙．

或者，可以在proActiv Github存储库中提出问题:https://github.com/GoekeLab/proActiv．

高伦雅芙为由

如果您使用proActiv，请注明:

demirciozylu, Deniz等，“泛癌症转录组分析揭示了通过替代启动子的普遍调控。”细胞178.6(2019): 1465-1477。

会话信息

#> R开发中(不稳定)(2022-10-25 r83175) #>平台:x86_64-pc-linux-gnu(64位)#>运行在Ubuntu 22.04.1 LTS下#> #>矩阵产品:默认#> BLAS: /home/biocbuild/bbs-3.17-bioc/R/lib/libRblas。so #> LAPACK: /usr/lib/x86_64-linux-gnu/ LAPACK /liblapack.so.3.10.0 #> #> locale: #> [1] LC_CTYPE=en_US。UTF-8 LC_NUMERIC= c# > [3] LC_TIME=en_GB LC_COLLATE= c# > [5] LC_MONETARY=en_US。utf - 8 LC_MESSAGES = en_US。UTF-8 #> [7] LC_PAPER=en_US。UTF-8 LC_NAME= c# > [9] LC_ADDRESS=C LC_TELEPHONE= c# > [11] LC_MEASUREMENT=en_US。UTF-8 LC_IDENTIFICATION= c# b> # b>附加基础包:#> [1]stats graphics grDevices utils datasets methods base #> #>其他附加包:#> [1]Rtsne_0.16 ggplot2_3.3.6 gridExtra_2.3 proActiv_1.9.0 #> #>通过命名空间加载(且未附加):# > [1] DBI_1.1.3 bitops_1.0-7 # > [3] biomaRt_2.55.0 rlang_1.0.6 # > [5] magrittr_2.0.3 matrixStats_0.62.0 # > [7] compiler_4.3.0 RSQLite_2.2.18 # > [9] GenomicFeatures_1.51.0 png_0.1-7 # > [11] vctrs_0.5.0 stringr_1.4.1 # > [13] pkgconfig_2.0.3 crayon_1.5.2 # > [15] fastmap_1.1.0 dbplyr_2.2.1 # > [17] XVector_0.39.0 ellipsis_0.3.2 # > [19] labeling_0.4.2 caTools_1.18.2 # > [21] utf8_1.2.2 Rsamtools_2.15.0 # > [23] rmarkdown_2.17 bit_4.0.4 # > [25] xfun_0.34 zlibbioc_1.45.0 # > [27] cachem_1.0.6GenomeInfoDb_1.35.0 # > [29] jsonlite_1.8.3 progress_1.2.2 # > [31] blob_1.2.3 highr_0.9 # > [33] DelayedArray_0.25.0 BiocParallel_1.33.0 # > [35] parallel_4.3.0 prettyunits_1.1.1 # > [37] R6_2.5.1 bslib_0.4.0 # > [39] stringi_1.7.8 RColorBrewer_1.1-3 # > [41] rtracklayer_1.59.0 genefilter_1.81.0 # > [43] GenomicRanges_1.51.0 jquerylib_0.1.4 # > [45] Rcpp_1.0.9 assertthat_0.2.1 # > [47] SummarizedExperiment_1.29.0 knitr_1.40 # > [49] IRanges_2.33.0 Matrix_1.5-1 # > [51] splines_4.3.0 tidyselect_1.2.0 # >[53] yaml_2.3.6 gplots_3.1.3 # > [55] codetools_0.2-18 curl_4.3.3 # > [57] lattice_0.20-45 tibble_3.1.8 # > [59] withr_2.5.0 Biobase_2.59.0 # > [61] KEGGREST_1.39.0 evaluate_0.17 # > [63] survival_3.4-0 BiocFileCache_2.7.0 # > [65] xml2_1.3.3 Biostrings_2.67.0 # > [67] pillar_1.8.1 filelock_1.0.2 # > [69] MatrixGenerics_1.11.0 KernSmooth_2.23-20 # > [71] stats4_4.3.0 generics_0.1.3 # > [73] rcurl_1.98 - 1.9 S4Vectors_0.37.0 # > [75] hms_1.1.2 munsell_0.5.0 # > [77] scales_1.2.1 gtools_3.9.3 # > [79]xtable_1.8-4 glue_1.6.2 #> [81] tools_4.3.0 bioplot_1 .9.0 #> [85] locfit_1.5-9.6 genome alignments_1.35.0 #> [87] XML_3.99-0.12 grid_4.3.0 #> [89] AnnotationDbi_1.61.0 colorspace_2.0-3 #> [91] GenomeInfoDbData_1.2.9 restfulr_0.0.15 #> [93] cli_3.4.1 rappdirs_0.3.3 #> [95] fansi_1.0.3 dplyr_1.0.10 #> [97] gtable_0.3.1 DESeq2_1.39.0 #> [99] sass_0.4.2 digest_0.6.30 #> [101] BiocGenerics_0.45.0 farver_2.1.1 #> [103] rjson_0.2.21 geneplotter_1.77.0 #> [105]Memoise_2.0.1 htmltools_0.5.3 #> [107] lifecycle_1.0.3 httr_1.4.4 #> [109] bit64_4.0.5