forked from duty-machine/duty-machine
-
Notifications
You must be signed in to change notification settings - Fork 32
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
- Loading branch information
Showing
1 changed file
with
15 additions
and
0 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,15 @@ | ||
--- | ||
title: "如何获取基因的转录起始/终止位点(TSS/TES)?" | ||
date: 2023-10-22T09:19:01Z | ||
draft: ["false"] | ||
tags: [ | ||
"fetched", | ||
"老俊俊的生信笔记" | ||
] | ||
categories: ["Acdemic"] | ||
--- | ||
如何获取基因的转录起始/终止位点(TSS/TES)? by 老俊俊的生信笔记 | ||
------ | ||
<div><section data-tool="mdnice编辑器" data-website="https://www.mdnice.com" data-mpa-powered-by="yiban.io"><h4 data-tool="mdnice编辑器"><span></span></h4><section><mp-common-profile data-pluginname="mpprofile" data-id="MzkyMTI1MTYxNA==" data-headimg="http://mmbiz.qpic.cn/mmbiz_png/G5jjcE4usey42oX5qyLTVibLRO9dz8ic5G4TpEHQc9rICYlpS4MHg6Et8cgXrQDqibvibXombicTro8t9cekJRlDBcw/0?wx_fmt=png" data-nickname="老俊俊的生信笔记" data-alias="JunJunLab" data-signature="老俊俊的生信技能和知识分享,我不是巨人,但你可以站在我的肩膀上更进一步!" data-from="0" data-is_biz_ban="0"></mp-common-profile></section><h4 data-tool="mdnice编辑器"><span></span></h4><h2 data-tool="mdnice编辑器"><span><span>1</span></span><span>引言</span><span></span></h2><blockquote data-tool="mdnice编辑器"><p>群里看到有粉丝卡在了怎么获取基因的转录起始位点 <strong>TSS</strong>数据上,导致后续画图发生错误。下面介绍一下如何获取基因的转录起始位点 <strong>TSS</strong>数据 以及转录终止位点 <strong>TES</strong> 数据。</p></blockquote><h2 data-tool="mdnice编辑器"><span><span>2</span></span><span>示例</span><span></span></h2><h3 data-tool="mdnice编辑器"><span></span><span>从 <strong>TXDB</strong> 注释文件获取</span><span></span></h3><p data-tool="mdnice编辑器"><strong>Bioconductor</strong> 上面有很多模式物种的注释包,你可以下载直接使用:</p><blockquote data-tool="mdnice编辑器"><p><strong>https://www.bioconductor.org/packages/release/BiocViews.html#___TxDb</strong></p></blockquote><figure data-tool="mdnice编辑器"><img data-ratio="0.5722108145106092" data-src="https://mmbiz.qpic.cn/sz_mmbiz_png/G5jjcE4usexmLWZ8rn8G1x7mp8MMyTIAId0zbRmicyZuicnJJ6z9AibUMg0Gqx1CQ9nAY3uVFpibfZ9FicbZic4dFalw/640?wx_fmt=png" data-type="png" data-w="1461" src="https://mmbiz.qpic.cn/sz_mmbiz_png/G5jjcE4usexmLWZ8rn8G1x7mp8MMyTIAId0zbRmicyZuicnJJ6z9AibUMg0Gqx1CQ9nAY3uVFpibfZ9FicbZic4dFalw/640?wx_fmt=png"></figure><p data-tool="mdnice编辑器">我们用小鼠的测试,首先获取基因位置信息:</p><blockquote data-tool="mdnice编辑器"><p>如果你要画整个基因上信号富集的信号的话,下面这个就可以直接给 <strong>normalizedMatrix/getTagMatrix</strong> 作为输入了。<strong>此外要注意的是你的 bigwig 文件的染色体前缀 chr 和你的 TSS 文件的 chr 保持一致,要么都有,要么都没有。不然后面做的也还是错的,一定要仔细!</strong></p></blockquote><pre data-tool="mdnice编辑器"><span></span><code><span>library</span>(tidyverse)<br><span>library</span>(rtracklayer)<br><span>library</span>(IRanges)<br><span>library</span>(TxDb.Mmusculus.UCSC.mm10.knownGene)<br><br>gene <- genes(TxDb.Mmusculus.UCSC.mm10.knownGene)<br>gene<br><span># GRanges object with 24528 ranges and 1 metadata column:</span><br><span># seqnames ranges strand | gene_id</span><br><span># <Rle> <IRanges> <Rle> | <character></span><br><span># 100009600 chr9 21062393-21073096 - | 100009600</span><br><span># 100009609 chr7 84935565-84964115 - | 100009609</span><br><span># 100009614 chr10 77711457-77712009 + | 100009614</span><br><span># ...</span><br><span># -------</span><br><span># seqinfo: 66 sequences (1 circular) from mm10 genome</span><br></code></pre><blockquote data-tool="mdnice编辑器"><p>基因相对于双链 DNA 是区分正负链的,所以 <strong>TSS/TES</strong> 位点就是基因位置的两边位置,只不过正负链基因相反而已。举个例子: 上面 <strong>100009600</strong> 这个基因在负链上,所以 <strong>TSS/TES</strong> 分别是 <strong>21073096/21062393</strong>。而 <strong>100009614</strong> 基因在正链上, 所以 <strong>TSS/TES</strong> 分别是 <strong>77711457/77712009</strong>。</p></blockquote><p data-tool="mdnice编辑器">我们使用 <strong>promoters</strong> 函数可以获取 TSS:</p><pre data-tool="mdnice编辑器"><span></span><code>tss <- promoters(gene,upstream = <span>0</span>,downstream = <span>1</span>)<br>tss<br><span># GRanges object with 24528 ranges and 1 metadata column:</span><br><span># seqnames ranges strand | gene_id</span><br><span># <Rle> <IRanges> <Rle> | <character></span><br><span># 100009600 chr9 21073096 - | 100009600</span><br><span># 100009609 chr7 84964115 - | 100009609</span><br><span># 100009614 chr10 77711457 + | 100009614</span><br><span># ...</span><br></code></pre><p data-tool="mdnice编辑器">我们也可以手动来算一下,很简单,思路就是正链取 <strong>start</strong>,负链取 <strong>end</strong> 即可, 可以看到结果和上面一致:</p><pre data-tool="mdnice编辑器"><span></span><code>gene_df <- data.frame(gene)<br><br><span># check</span><br>head(gene_df,<span>3</span>)<br><span># seqnames start end width strand gene_id</span><br><span># 1 chr9 21062393 21073096 10704 - 100009600</span><br><span># 2 chr7 84935565 84964115 28551 - 100009609</span><br><span># 3 chr10 77711457 77712009 553 + 100009614</span><br><br>tss2 <- gene_df |><br> mutate(start = ifelse(strand == <span>"+"</span>,start,end),<br> end = start) |><br> GRanges()<br>tss2<br><span># GRanges object with 24528 ranges and 1 metadata column:</span><br><span># seqnames ranges strand | gene_id</span><br><span># <Rle> <IRanges> <Rle> | <character></span><br><span># 100009600 chr9 21073096 - | 100009600</span><br><span># 100009609 chr7 84964115 - | 100009609</span><br><span># 100009614 chr10 77711457 + | 100009614</span><br><span># ...</span><br></code></pre><hr data-tool="mdnice编辑器"><p data-tool="mdnice编辑器">思路反过来,我们就可以获取 <strong>TES</strong> 的信息了:</p><pre data-tool="mdnice编辑器"><span></span><code>tes <- gene_df |><br> mutate(start = ifelse(strand == <span>"+"</span>,end,start),<br> end = start) |><br> GRanges()<br>tes<br><span># GRanges object with 24528 ranges and 1 metadata column:</span><br><span># seqnames ranges strand | gene_id</span><br><span># <Rle> <IRanges> <Rle> | <character></span><br><span># [1] chr9 21062393 - | 100009600</span><br><span># [2] chr7 84935565 - | 100009609</span><br><span># [3] chr10 77712009 + | 100009614</span><br><span># ...</span><br></code></pre><h3 data-tool="mdnice编辑器"><span></span><span>从 GTF 文件获取</span><span></span></h3><p data-tool="mdnice编辑器">有时候可能没有你的物种注释好的包,你可能只有注释 <strong>GTF</strong> 注释文件,根据上面的思路,我们也可以从 <strong>GTF</strong> 文件来准备。</p><pre data-tool="mdnice编辑器"><span></span><code>gtf <- import.gff(<span>"../Mus_musculus.GRCm38.102.gtf.gz"</span>,format = <span>"gtf"</span>) |><br> data.frame()<br></code></pre><figure data-tool="mdnice编辑器"><img data-ratio="0.2681734880879658" data-src="https://mmbiz.qpic.cn/sz_mmbiz_png/G5jjcE4usexmLWZ8rn8G1x7mp8MMyTIABFheibibuUUy8H7ibNc1vPNaB50clTdk8Riasbd0oichhEnNczaibIqVficog/640?wx_fmt=png" data-type="png" data-w="1637" src="https://mmbiz.qpic.cn/sz_mmbiz_png/G5jjcE4usexmLWZ8rn8G1x7mp8MMyTIABFheibibuUUy8H7ibNc1vPNaB50clTdk8Riasbd0oichhEnNczaibIqVficog/640?wx_fmt=png"></figure><p data-tool="mdnice编辑器">信息很多,提取我们需要的即可:</p><pre data-tool="mdnice编辑器"><span></span><code>df_gene <- gtf |><br> dplyr::filter(type == <span>"gene"</span>) |><br> dplyr::select(seqnames,start,end,gene_id,gene_name,strand)<br><br><span># check</span><br>head(df_gene,<span>3</span>)<br><span># seqnames start end gene_id gene_name strand</span><br><span># 1 1 3073253 3074322 ENSMUSG00000102693 4933401J01Rik +</span><br><span># 2 1 3102016 3102125 ENSMUSG00000064842 Gm26206 +</span><br><span># 3 1 3205901 3671498 ENSMUSG00000051951 Xkr4 -</span><br></code></pre><p data-tool="mdnice编辑器">到这里了方法就和上面提取的一样了,这里就不再写了。这里你可以看到染色体名称是没有 <strong>chr</strong> 前缀的,你需要检查和你的 <strong>bigwig</strong> 文件是否保持一致,才好做后面的分析!</p><h2 data-tool="mdnice编辑器"><span><span>3</span></span><span>嫌麻烦? 打包成函数</span><span></span></h2><blockquote data-tool="mdnice编辑器"><p>如果你还嫌麻烦,那就把上面代码打包成一个小函数,这样方便快捷。</p></blockquote><p data-tool="mdnice编辑器">示例:</p><pre data-tool="mdnice编辑器"><span></span><code>getSite(object = TxDb.Mmusculus.UCSC.mm10.knownGene,type = <span>"tss"</span>)<br><br>getSite(object = TxDb.Mmusculus.UCSC.mm10.knownGene,type = <span>"tss"</span>,<br> upstream = <span>10</span>,downstream = <span>5</span>)<br> <br>getSite(object = TxDb.Mmusculus.UCSC.mm10.knownGene,type = <span>"tes"</span>)<br><br>getSite(object = TxDb.Mmusculus.UCSC.mm10.knownGene,type = <span>"body"</span>)<br><br>getSite(object = <span>"../Mus_musculus.GRCm38.102.gtf.gz"</span>,type = <span>"tss"</span>)<br></code></pre><p data-tool="mdnice编辑器">函数代码:</p><pre data-tool="mdnice编辑器"><span></span><code>getSite <- <span>function</span>(object = <span>NULL</span>,<br> upstream = <span>0</span>,<br> downstream = <span>0</span>,<br> type = c(<span>"tss"</span>,<span>"tes"</span>,<span>"body"</span>)){<br> type <- match.arg(type,c(<span>"tss"</span>,<span>"tes"</span>,<span>"body"</span>))<br><br> <span># ============================================================================</span><br> <span># define func</span><br> <span># ============================================================================</span><br> getfun <- <span>function</span>(data = <span>NULL</span>,type = <span>NULL</span>,upstream = <span>0</span>,downstream = <span>0</span>){<br> <span># check type to extract</span><br> <span>if</span>(type == <span>"tss"</span>){<br> ts <- data |><br> dplyr::mutate(start = ifelse(strand == <span>"+"</span>,start,end),end = start)<br><br> }<span>else</span> <span>if</span>(type == <span>"tes"</span>){<br> ts <- data |><br> dplyr::mutate(start = ifelse(strand == <span>"+"</span>,end,start),end = start)<br> }<span>else</span>{<br> ts <- data<br> }<br><br> <span># extend</span><br> ts_ed <- ts |><br> dplyr::mutate(start1 = ifelse(strand == <span>"+"</span>,start - upstream,start - downstream),<br> end1 = ifelse(strand == <span>"+"</span>,end + downstream,end + upstream)) |><br> dplyr::select(-start,-end) |><br> dplyr::rename(start = start1,end = end1) |><br> GRanges()<br><br> <span>return</span>(ts_ed)<br> }<br><br> <span># ============================================================================</span><br> <span># extract from file</span><br> <span># ============================================================================</span><br> <span>if</span>(<span>"TxDb"</span> %<span>in</span>% class(object)){<br> gene <- GenomicFeatures::genes(object)<br> gene_df <- data.frame(gene)<br><br> res <- getfun(data = gene_df,type = type,<br> upstream = upstream,downstream = downstream)<br><br> }<span>else</span> <span>if</span>(is.character(object)){<br> gtf <-rtracklayer::import.gff(object,format = <span>"gtf"</span>) |><br> data.frame()<br><br> ctype <- ifelse(<span>"gene"</span> %<span>in</span>% unique(gtf$type),<span>"gene"</span>,<span>"transcript"</span>)<br><br> gene_df <- gtf |><br> dplyr::filter(type == ctype) |><br> dplyr::select(seqnames,start,end,gene_id,gene_name,strand)<br><br> res <- getfun(data = gene_df,type = type,<br> upstream = upstream,downstream = downstream)<br> }<br><br> <span>return</span>(res)<br>}<br></code></pre><h2 data-tool="mdnice编辑器"><span><span>4</span></span><span>结尾</span><span></span></h2><blockquote data-tool="mdnice编辑器"><p><strong>路漫漫其修远兮,吾将上下而求索。</strong></p></blockquote><hr data-tool="mdnice编辑器"><p data-tool="mdnice编辑器">欢迎加入生信交流群。加我微信我也拉你进 <strong>微信群聊</strong> <strong>老俊俊生信交流群</strong> <strong>(微信交流群需收取 20 元入群费用,一旦交费,拒不退还!(防止骗子和便于管理))</strong> 。QQ 群可免费加入, 记得进群按格式修改备注哦。</p><section data-tool="mdnice编辑器"><section><p><strong>老俊俊微信:</strong></p><figure><img data-ratio="1" data-src="https://mmbiz.qpic.cn/sz_mmbiz_png/G5jjcE4usexmLWZ8rn8G1x7mp8MMyTIAMAjv4WFN2yRnVzHa5Aq4xl3bRce9aXyvsC3giaL6iaHg5eI18gJsjSWg/640?wx_fmt=png" data-type="png" data-w="430" src="https://mmbiz.qpic.cn/sz_mmbiz_png/G5jjcE4usexmLWZ8rn8G1x7mp8MMyTIAMAjv4WFN2yRnVzHa5Aq4xl3bRce9aXyvsC3giaL6iaHg5eI18gJsjSWg/640?wx_fmt=png"></figure><figure><img data-ratio="1.3668430335097002" data-src="https://mmbiz.qpic.cn/sz_mmbiz_png/G5jjcE4usexmLWZ8rn8G1x7mp8MMyTIAbM7l2cmGQkXgMtBLLU5ODHWzBo0pgfoOMhkmzStHONIukPbxgEyLSA/640?wx_fmt=png" data-type="png" data-w="567" src="https://mmbiz.qpic.cn/sz_mmbiz_png/G5jjcE4usexmLWZ8rn8G1x7mp8MMyTIAbM7l2cmGQkXgMtBLLU5ODHWzBo0pgfoOMhkmzStHONIukPbxgEyLSA/640?wx_fmt=png"></figure></section><section><p><strong>知识星球:</strong></p><figure><img data-ratio="1.5896226415094339" data-src="https://mmbiz.qpic.cn/sz_mmbiz_jpg/G5jjcE4usexmLWZ8rn8G1x7mp8MMyTIANerIBOfCR2fMuicyY73ZcD5QLUW1J6ticzdlWsZvv5wM3c5MroIHu6sA/640?wx_fmt=jpeg" data-type="jpeg" data-w="1060" src="https://mmbiz.qpic.cn/sz_mmbiz_jpg/G5jjcE4usexmLWZ8rn8G1x7mp8MMyTIANerIBOfCR2fMuicyY73ZcD5QLUW1J6ticzdlWsZvv5wM3c5MroIHu6sA/640?wx_fmt=jpeg"></figure></section></section><hr data-tool="mdnice编辑器"><p data-tool="mdnice编辑器"><br></p><center data-tool="mdnice编辑器"><strong> 往期回顾目录</strong></center><blockquote data-tool="mdnice编辑器"><p><br></p><center><strong><a href="https://mp.weixin.qq.com/s?__biz=MzkyMTI1MTYxNA==&mid=2247510151&idx=1&sn=de70b667d6eb7f581a68fd3fc231b1fc&chksm=c18492f6f6f31be05a57cb681d802fb813e36ac67f3d70fb05d26a6ea1aeee14bd17a91844e7&token=1624839730&lang=zh_CN&scene=21#wechat_redirect" data-linktype="2">ChIP-seq 绘图设计</a></strong></center><strong><center><a href="https://mp.weixin.qq.com/s?__biz=MzkyMTI1MTYxNA==&mid=2247510122&idx=1&sn=f06cfbecb15ee4f133c9b3d49791743a&chksm=c184921bf6f31b0d471e37ef63fdc4c7caf63df2a6598bd1f3c04d9b1505e2e103ef4f60eb66&token=1598639854&lang=zh_CN&scene=21#wechat_redirect" data-linktype="2">ggChIPvis 可视化基因组富集信号</a></center></strong><strong><center><a href="https://mp.weixin.qq.com/s?__biz=MzkyMTI1MTYxNA==&mid=2247510122&idx=1&sn=f06cfbecb15ee4f133c9b3d49791743a&chksm=c184921bf6f31b0d471e37ef63fdc4c7caf63df2a6598bd1f3c04d9b1505e2e103ef4f60eb66&token=1598639854&lang=zh_CN&scene=21#wechat_redirect" data-linktype="2">ggChIPvis 可视化基因组富集信号</a></center></strong><strong><center><a href="https://mp.weixin.qq.com/s?__biz=MzkyMTI1MTYxNA==&mid=2247510110&idx=1&sn=5b8a12b864509d65c243b61d080350a9&chksm=c184922ff6f31b3943a3681def268998388a3a930b2c384f3d8ea33ff17a437cc4513fc1b286&token=819162780&lang=zh_CN&scene=21#wechat_redirect" data-linktype="2">ChIP-seq 可视化探索</a></center></strong><strong><center><a href="https://mp.weixin.qq.com/s?__biz=MzkyMTI1MTYxNA==&mid=2247510080&idx=1&sn=d296de1de6d26ca85c724c392fe9959b&chksm=c1849231f6f31b272f3826f9a8f5277a1272f44443a5553513c0ffe3d82ca48d3c518ad49d0f&token=1127717108&lang=zh_CN&poc_token=HItqJWWjkYgH0feCwRqHpggUnZg-boP69onpft9F&scene=21#wechat_redirect" data-linktype="2">自定义基因集 GseaVis 绘图</a></center></strong><strong><center><a href="https://mp.weixin.qq.com/s?__biz=MzkyMTI1MTYxNA==&mid=2247510066&idx=1&sn=19c713554cd00a6d444e098145f3d286&chksm=c1849243f6f31b5582896d837ce91f6655d7029fc6c95bfc98c622787d4a9dfe57a74152fd75&token=1127717108&lang=zh_CN&scene=21#wechat_redirect" data-linktype="2">批量下载 starbase 数据</a></center></strong><strong><center><a href="https://mp.weixin.qq.com/s?__biz=MzkyMTI1MTYxNA==&mid=2247510050&idx=1&sn=9b1a5fb1cc6ecf52844847f4aba6bfb8&chksm=c1849253f6f31b450c12a6748a188e734f92c16be1abc51d75afabb33f404f45d60d1d616191&token=100163832&lang=zh_CN&scene=21#wechat_redirect" data-linktype="2">听说你想在 python 中对 kegg 可视化?</a></center></strong><strong><center><a href="https://mp.weixin.qq.com/s?__biz=MzkyMTI1MTYxNA==&mid=2247510034&idx=1&sn=e43a22b9de87e7e4500d3e2f52bc50c5&chksm=c1849263f6f31b75338900b88e5abfe2747e7cd980706c8348e23078d65ca28273cbcbf2ff84&token=203556469&lang=zh_CN&scene=21#wechat_redirect" data-linktype="2">优雅的使用 lattice 可视化单细胞</a></center></strong><strong><center><a href="https://mp.weixin.qq.com/s?__biz=MzkyMTI1MTYxNA==&mid=2247510008&idx=1&sn=86892dd2412183414ed984d5279848fc&chksm=c1849d89f6f3149f5771e1b74a7e9ad1e72c619c203e227f25bd895712c981c2f9d8230e91b6&token=203556469&lang=zh_CN&scene=21#wechat_redirect" data-linktype="2">听说你想在 circlize 的环形热图里添加显著性?</a></center></strong><strong><center><a href="https://mp.weixin.qq.com/s?__biz=MzkyMTI1MTYxNA==&mid=2247509978&idx=1&sn=7af95b9220830dfc92a31ef64092e72e&chksm=c1849dabf6f314bd5d28fee0a1467696de2c157d506ea1d19a43c2319811514fcc76515469bf&token=1162368673&lang=zh_CN&scene=21#wechat_redirect" data-linktype="2">听说表达趋势你想分面展示?</a></center></strong><p><br></p></blockquote></section><p><mp-style-type data-value="3"></mp-style-type></p></div> | ||
<hr> | ||
<a href="https://mp.weixin.qq.com/s/cUbEDt6P46YFPI1OkyuuHw",target="_blank" rel="noopener noreferrer">原文链接</a> |