当前位置: 首页 > news >正文

百度自然排名网站的logo怎么做长春做网站推广的公司

百度自然排名网站的logo怎么做,长春做网站推广的公司,网站建设员招聘,皮具 东莞网站建设在正式分析之前,对于数据的处理是至关重要的,这种重要性是体现在很多方面,其中有一点是要求分析者采用正确的数据类型。 对于芯片数据,原始数据进行log2处理之后可以进行很多常见的分析,比如差异分析、热图、箱线图、…

在正式分析之前,对于数据的处理是至关重要的,这种重要性是体现在很多方面,其中有一点是要求分析者采用正确的数据类型。

对于芯片数据,原始数据进行log2处理之后可以进行很多常见的分析,比如差异分析、热图、箱线图、PCA分析、生存分析、模型构建,聚类分析和相关性分析等。

对于转录组数据,在上述的常见分析中只有差异分析时需要采用count值其他的分析是需要采用log2后的cpm,tpm,fpkm,rpkm数据

除了上述的常规分析,在使用不同R包进行分析之前务必浏览一下输入数据的要求

那么芯片数据还好说,毕竟后续进行log2处理后就可以做很多分析。但是转录组数据的可选项就比较多了。

但目前常用的只有tpm和cpm

count数据转换为cpm数据非常简单

# exprSet是count表达矩阵
# 一句代码搞定
exprSet = log2(edgeR::cpm(exprSet)+1)

比较难的就是count数据转换为tpm数据,因此搬运了常规的流程和R包的方法,做个对比

首先要去获取基因长度文件,因为后续需要用这个数据去矫正基因长度。

网址:https://gdc.cancer.gov/about-data/gdc-data-processing/gdc-reference-files

方法一:常规方法-来自生信技能树

1、获取gencode.v36.annotation并处理

if(!require(GenomicFeatures))BiocManager::install("GenomicFeatures")
library(GenomicFeatures)
txdb <- makeTxDbFromGFF("gencode.v36.annotation.gtf.gz",format="gtf")
# 获取每个基因id的外显子数据
exons.list.per.gene <- exonsBy(txdb,by="gene")
# 对于每个基因,将所有外显子减少成一组非重叠外显子,计算它们的长度(宽度)并求和
exonic.gene.sizes <- sum(width(GenomicRanges::reduce(exons.list.per.gene)))
# 得到geneid和长度数据
gfe <- data.frame(gene_id=names(exonic.gene.sizes),length=exonic.gene.sizes)
head(gfe)[1:5,1:2]
#                               gene_id length
# ENSG00000000003.15 ENSG00000000003.15   4536
# ENSG00000000005.6   ENSG00000000005.6   1476
# ENSG00000000419.13 ENSG00000000419.13   1207
# ENSG00000000457.14 ENSG00000000457.14   6883
# ENSG00000000460.17 ENSG00000000460.17   5970
save(gfe,file = "gfe.Rdata")

2、使用TCGA-HNSC数据

load("hnsc_exp.Rdata")
head(hnsc)[1:3,1:3]
#                      TCGA-UF-A7JF-01A-11R-A34R-07 TCGA-CN-4725-01A-01R-1436-07 TCGA-D6-6827-01A-11R-1915-07
# ENSG00000000003.15                         2090                         1680                         5433
# ENSG00000000005.6                             0                            0                            0
# ENSG00000000419.13                         2098                         3872                         2240
identical(rownames(hnsc),rownames(gfe))
#[1] TRUE 行名是能够对上的

使用了TCGA-HNSC中的count数据,检测了一下count数据和下载的基因信息的顺序是一致的。

3、曾老师写的代码进行count/tpm转化

load("gfe.Rdata")
#提取基因长度列
effLen = gfe$length
#转化
Counts2TPM <- function(counts, effLen){rate <- log(counts) - log(effLen)denom <- log(sum(exp(rate)))exp(rate - denom + log(1e6))
}
hnsc_tpm_raw <- apply(hnsc, 2, Counts2TPM, effLen = effLen)
head(hnsc_tpm_raw)[1:3,1:3]
#                      TCGA-UF-A7JF-01A-11R-A34R-07 TCGA-CN-4725-01A-01R-1436-07 TCGA-D6-6827-01A-11R-1915-07
# ENSG00000000003.15                     31.13713                     15.74248                     72.37614
# ENSG00000000005.6                       0.00000                      0.00000                      0.00000
# ENSG00000000419.13                    117.46366                    136.35307                    112.14231

4、加载既往处理好的数据进行对比

load("hnsc_tpm.Rdata")
head(hnsc_tpm)[1:3,1:3]
#                       TCGA-UF-A7JF-01A-11R-A34R-07 TCGA-CN-4725-01A-01R-1436-07 TCGA-D6-6827-01A-11R-1915-07
# ENSG00000000003.15                     31.13713                     15.74248                     72.37614
# ENSG00000000005.6                       0.00000                      0.00000                      0.00000
# ENSG00000000419.13                    117.46366                    136.35307                    112.14231

结果是一致的

方法二:R包 DGEobj.utils
load("gfe.Rdata")
#提取上面处理好的基因长度列
effLen = gfe$length
library(DGEobj.utils)
CC_res <- convertCounts(exp,unit = "TPM", # CPM、FPKM、FPK 或 TPMgeneLength = effLen, #这里还是需要下在基因长度数据log = FALSE, #默认为False,设为TRUE时返回Log2值normalize = "none", #默认为none,调用edgeR的calcNormFactors()进行标准化prior.count = NULL #为避免取0的对数,对每个观测值添加平均count。仅当 log = TRUE 时使用,
)
head(CC_res)[1:3,1:3]
#                   TCGA-UF-A7JF-01A-11R-A34R-07 TCGA-CN-4725-01A-01R-1436-07 TCGA-D6-6827-01A-11R-1915-07
#ENSG00000000003.15                     31.13713                     15.74248                     72.37614
#ENSG00000000005.6                       0.00000                      0.00000                      0.00000
#ENSG00000000419.13                    117.46366                    136.35307                    112.14231

与前面得到的分析结果是一致的。

需要下载基因长度数据,但是前期处理完之后后面可以很方便的转化为各种想要的格式(CPM、FPKM、FPK 或 TPM)。

方法三:R包 IOBR
load("hnsc_exp.Rdata")
exp <- hnsc
library(IOBR)
exp_tpm <- count2tpm(exp, #表达矩阵, 行为基因idType = "Ensembl", #"Ensembl" "ENTREZ","SYMBOL" 
)
head(exp_tpm)[1:3,1:3]
#         TCGA-UF-A7JF-01A-11R-A34R-07 TCGA-CN-4725-01A-01R-1436-07 TCGA-D6-6827-01A-11R-1915-07
# TSPAN6                     31.89988                     16.16771                     73.85059
# TNMD                        0.00000                      0.00000                      0.00000
# DPM1                       15.65888                     18.22162                     14.88931
# 可以看到这个函数还自动把gene_id给转换成了symbol# 为了更好的对比,我们也把hnsc_tpm_raw中的gene_id转换成symbol
# 使用小洁老师的trans_exp_new函数
library(tinyarray)
a <- trans_exp_new(hnsc_tpm_raw)
head(a)[1:3,1:3]
#           TCGA-UF-A7JF-01A-11R-A34R-07 TCGA-CN-4725-01A-01R-1436-07 TCGA-D6-6827-01A-11R-1915-07
# DDX11L1                       0.000000                    0.0000000                     0.000000
# WASH7P                        1.150477                    0.8180029                     2.191641
# MIR6859-1                     0.993794                    0.0000000                     4.443138# 排个序再看一下
p <- identical(rownames(a),rownames(exp_tpm));p
if(!p) {s = intersect(rownames(a),rownames(exp_tpm))a = a[s,]exp_tpm = exp_tpm[s,]
}
head(a)[1:3,1:3]
#          TCGA-UF-A7JF-01A-11R-A34R-07 TCGA-CN-4725-01A-01R-1436-07 TCGA-D6-6827-01A-11R-1915-07
# DDX11L1                     0.000000                   0.00000000                     0.000000
# WASH7P                      1.150477                   0.81800290                     2.191641
# FAM138A                     0.000000                   0.03486849                     0.000000head(exp_tpm)[1:3,1:3]
#          TCGA-UF-A7JF-01A-11R-A34R-07 TCGA-CN-4725-01A-01R-1436-07 TCGA-D6-6827-01A-11R-1915-07
# DDX11L1                      0.00000                   0.00000000                      0.00000
# WASH7P                      22.63327                   5.95811827                     15.88102
# FAM138A                      0.00000                   0.03581036                      0.00000

IOBR包中的count2tpm函数只能进行tpm转化(github上搜了这个R包内容似乎没有转化为其他数据格式的函数了)。用默认的方式进行运算得到的结果存在一定的偏差,而且我个人觉得这个偏差有点大... 但是我暂时不知道是什么原因?是内置的基因长度顺序有问题?还是我某个参数设的不对?求高手指点~

综合上述分析,暂时还是选择常规方法DGEobj.utils R包中的convertCounts函数吧~

参考资料:

1、https://hbctraining.github.io/Training-modules/planning_successful_rnaseq/lessons/sample_level_QC.html

2、https://rdrr.io/cran/DGEobj.utils/man/convertCounts.html

3、https://rdrr.io/github/IOBR/IOBR/man/count2tpm.html

4、生信技能树推文:https://mp.weixin.qq.com/s/IUV9dSbRBK1nvetixKOCRw

致谢:感谢曾老师,小洁老师以及生信技能树团队全体成员。

:若对内容有疑惑或者有发现明确错误的朋友,请联系后台(欢迎交流)。更多内容可关注公众号:生信方舟

- END -

http://www.wangmingla.cn/news/21390.html

相关文章:

  • 重庆做网站公司哪家好郑州抖音推广
  • 上海公司注册网站打开百度网站
  • 网站开发保密合同范本百度指数分析数据
  • 网站建设租用服务器天津百度推广公司地址
  • 广州口碑好的网站建设定制免费网页在线客服系统代码
  • 百度广告联盟官网下载热狗网站关键词优化
  • 西部数码网站管理助手 卸载企业营销网站制作
  • 有园林案例的网站可以下载新闻视频的网站
  • 泉州网站建设技术外包谷歌关键词优化怎么做
  • java做的小说网站的源代码正规百度推广
  • 互联网公司网站建设ppt模板下载媒体广告投放平台
  • 宿迁做网站的百度seo官方网站
  • 百度关键词策划和seo的优化seo关键词排名优化系统
  • 合肥效果好的网站推广杭州排名推广
  • 郑州网站制作费用英文seo实战派
  • 定手机网站建设拉新项目官方一手平台
  • 做棋牌网站团队潮州网站建设
  • 摄影网站设计图片今天
  • 资兴做网站公司seo优化检测
  • 网站域名的管理密码如何索取网推怎么做最有效
  • 网站开发项目时序图购物网站网页设计
  • 湘潭建设厅官方网站百度识图网页版
  • 个人网站设计 优帮云今日新闻最新头条
  • 哪里可做网站石家庄网络seo推广
  • 定制企业网站怎样设计一个网页
  • sunjava做网站如何网站seo
  • 好看的美食怎么做视频网站比较好用的搜索引擎
  • 免费域名分发网站别做网络推广员
  • 智慧房产信息管理平台免费seo推广计划
  • 网站建设通知书网站结构优化