内容

安装包

安装所有必需的软件包可以通过以下命令完成:

library(BiocInstaller) biocLite(c("waldronlab/MicrobiomeWorkshop"， "epiviz/metavizr@bioc-workshop")， dependencies=TRUE)

仔细检查“metavizr”是从GitHub安装的:

如果(!("metavizr" %in% installed.packages()) || "1.1.3" != installed.packages()["metavizr"，"Version"]) devtools::install_github("epiviz/metavizr@bioc-workshop"， build_vignettes = TRUE)
##从URL https://api.github.com/repos/epiviz/metavizr/zipball/bioc-workshop下载GitHub repo epiviz/metavizr@bioc-workshop ##
##安装metavizr
## '/usr/local/lib/R/bin/R'——no-site-file——no-environ——no-save \ ##——no-restore——quiet CMD build \ ## '/tmp/RtmpF7ZUtk/devtools1d676e54f4b/epiviz-metavizr-f4884e2' \ ##——no-save -data——no-manual
# #
##安装失败:Command failed (1)

加载所有需要的软件包:

cranpkgs=c("ggplot2"，"devtools"，"vegan"，"httr") BioCpkgs=c("curatedMetagenomicData"， "phyloseq") ghpkgs= "metavizr" allpkg <- c(cranpkgs, BioCpkgs, ghpkgs) suppressPackageStartupMessages(sapply(allpkg, require, character.)only = TRUE))
在“RSQLite”中没有找到请求:dbGetQuery的方法
## ggplot2 devtools vegan ## TRUE TRUE ## httr curatedMetagenomicData phyloseq ## TRUE TRUE ## metavizr ## TRUE

curatedMetagenomicData

curatedMetagenomicData为每个数据集提供6种数据类型:

  1. 种水平的分类概况，用从界到品系水平的相对丰度表示
  2. 存在独特的进化支特异性标记
  3. 丰富的独特的，进化枝特异性标记
  4. 基因家族的丰度
  5. 代谢途径覆盖率
  6. 代谢途径丰度

类型1-3由MetaPhlAn2;4-6是由HUMAnN2．

目前,curatedMetagenomicData提供:

这些数据集记录在参考手册释放或重击版本的Bioconductor。开发版本有大约两倍的数据集和更多的功能。Bioc2017研讨会提供的AMI运行Bioconductor的开发版本;看到//www.andersvercelli.com/2021欧洲杯体育投注开户developers/how-to/useDevel/有关将您自己的机器升级到开发版本以及在开发和发布之间切换的说明。

使用curatedMetagenomicData资源

资源的使用curatedMetagenomicData通过使用Bioconductor的ExperimentHub平台简化了该过程，该平台允许通过直观的界面访问数据。首先,curatedMetagenomicData安装时使用BiocInstaller然后作为库调用——这个过程允许用户简单地将数据集作为函数调用，因为包知道存在于ExperimentHub S3桶中的资源。

# BiocInstaller::biocLite("curatedMetagenomicData") # BiocInstaller::biocLite("waldronlab/curatedMetagenomicData") #bleeding edge version suppressPackageStartupMessages(library(curatedMetagenomicData))

可用的示例和元数据

所有可用样本的手动管理元数据都在单个表中提供combined_metadata：

? combined_metadata视图(combined_metadata)
表(combined_metadata antibiotics_current_use美元)
## ##否是## 1855 558
表(combined_metadata疾病美元)
AR # # # # # #广告广告;10 16 # #广告;哮喘广告;哮喘;AR # 8 # 4 # # AR CDI; NA # 14 # 2 # # CDI;蜂窝织炎CDI;骨关节炎# 1 # 2 # # CDI;肺炎CDI; ureteralstone 15 # # 1 # # CRC CRC; T2D;高血压231 # # 2 # # CRC; fatty_liver CRC; fatty_liver; T2D;高血压3 # # 4 # # CRC; fatty_liver;高血压CRC;高血压12 # # 12 # # IBD IGT # # 148 49 # # NK STEC 43 # # 1 # #日内转内转;腹腔;IBD 24 # # 3 # # T2D T2D;高血压223 # # 1 # # abdominalhernia腺瘤# # 2 50 # #腺瘤;T2D腺瘤;T2D;高血压# #1 2 ## adenoma;fatty_liver adenoma;fatty_liver;T2D ## 12 1 ## adenoma;fatty_liver;T2D;hypertension adenoma;fatty_liver;hypertension ## 1 14 ## adenoma;hypertension arthritis ## 8 11 ## asthma;AR bronchitis ## 8 18 ## cellulitis cirrhosis ## 6 8 ## cirrhosis;HBV cirrhosis;HBV;HDV ## 39 3 ## cirrhosis;HBV;HDV;ascites cirrhosis;HBV;HE ## 2 2 ## cirrhosis;HBV;HEV cirrhosis;HBV;HEV;ascites ## 3 4 ## cirrhosis;HBV;ascites cirrhosis;HBV;schistosoma;ascites ## 44 1 ## cirrhosis;HBV;wilson;ascites cirrhosis;HCV ## 1 1 ## cirrhosis;HEV;HE cirrhosis;ascites ## 1 9 ## cirrhosis;pbcirrhosis;ascites cirrhosis;schistosoma;HEV;ascites ## 1 2 ## cirrhosis;schistosoma;ascites cirrhosis;wilson;ascites ## 1 1 ## coeliac;gestational_diabetes;CMV cough ## 2 2 ## cystitis fatty_liver ## 1 8 ## fatty_liver;T2D fatty_liver;T2D;hypertension ## 3 9 ## fatty_liver;hypertension fever ## 13 3 ## gangrene healthy ## 1 2455 ## hepatitis hypertension ## 3 11 ## infectiousgastroenteritis osteoarthritis ## 5 1 ## otitis pneumonia ## 107 7 ## psoriasis psoriasis;arthritis ## 36 12 ## pyelonefritis pyelonephritis ## 2 6 ## respiratoryinf salmonellosis ## 13 1 ## schizophrenia schizophrenia;CD ## 12 1 ## schizophrenia;T2D sepsis ## 3 1 ## skininf stomatitis ## 2 2 ## suspinf tonsillitis ## 1 3

访问数据集

单个数据项目可以通过每个数据集函数或curatedMetagenomicData ()函数。对数据集名称的函数调用返回一个Bioconductor ExpressionSet对象:

suppressPackageStartupMessages(library(curatedMetagenomicData)) zeller.eset = ZellerG_2014.metaphlan_bugs_list.stool()

下面的代码创建了一个包含两个元素的列表ExpressionSet提供BritoIL_2016和Castro-NallarE_2015口腔分类资料的对象(此方法适用于release或devel，但这些数据集仅适用于devel)。请参阅PDF手册中可用数据集的完整列表。

“BritoIL_2016.metaphlan_bugs_list.”或alcavity", "Castro-NallarE_2015.metaphlan_bugs_list.oralcavity") esl <- curatedMetagenomicData(oral, dryrun = FALSE) esl esl[[1]] esl[[2]]

下面将提供所有粪便数据集，如果dryrun = FALSE设置:

curatedMetagenomicData(“* metaphlan_bugs_list。大便*"，dryrun = TRUE)
##试运行:查看要下载的数据集的返回值。运行' dryrun=FALSE '来实际下载这些数据集。
##[1]”AsnicarF_2017.metaphlan_bugs_list。凳子“##[2]”BritoIL_2016.metaphlan_bugs_list。凳子" ## [3]"FengQ_2015.metaphlan_bugs_list。stool" ## [4] "HMP_2012.metaphlan_bugs_list。凳子" ## [5]"Heitz-BuschartA_2016.metaphlan_bugs_list。凳子“##[6]”KarlssonFH_2013.metaphlan_bugs_list。stool" ## [7] "LeChatelierE_2013.metaphlan_bugs_list。stool" ## [8] "LiuW_2016.metaphlan_bugs_list。stool" ## [9] "LomanNJ_2013.metaphlan_bugs_list。凳子“##[10]”尼尔森hb_2014 .metaphlan_bugs_list。stool" ## [11] "Obregon-TitoAJ_2015.metaphlan_bugs_list。凳子" ## [12]"QinJ_2012.metaphlan_bugs_list。凳子“##[13]”QinN_2014.metaphlan_bugs_list。stool" ## [14] "RampelliS_2015.metaphlan_bugs_list。stool" ## [15] "RaymondF_2016.metaphlan_bugs_list。凳子" ## [16]"SchirmerM_2016.metaphlan_bugs_list。凳子“##[17]”VatanenT_2016.metaphlan_bugs_list。凳子“##[18]”VincentC_2016.metaphlan_bugs_list。“##[19]”VogtmannE_2016.metaphlan_bugs_list。凳子“##[20]”谢h_2016 .metaphlan_bugs_list。凳子" ## [21]"YuJ_2015.metaphlan_bugs_list。stool" ## [22] "ZellerG_2014.metaphlan_bugs_list.stool"

合并多个数据集

下面将上面下载的两个口腔数据集合并为一个ExpressionSet(仅限devel)。

eset <- mergeData(esl) esset
## ExpressionSet (storageMode: lockkedenvironment) ## assayData: 1914个特征，172个样本##元素名称:exps ## protocolData: none ## phenoData ## sampleNames: BritoIL_2016.metaphlan_bugs_list.oralcavity:M1.1。## britoil_2016 .metaphlan_bugs_list.oral - cavity:M1.10。山……# # Castro-NallarE_2015.metaphlan_bugs_list。或alcavity:ES_080 (172 ## total) ## varLabels: body_site antibiotics_current_use ... studyID (19 ## total) ## varMetadata: labelDescription ## featureData: none ## experimentData: use 'experimentData(object)' ## Annotation:

这适用于任何数量的数据集。该函数不会阻止您合并不同的数据类型(例如，带有基因家族的翻译错误列表)，但您可能不想这样做。

使用ExpressionSet对象

所有数据集表示为ExpressionSet对象，因为类的集成特性及其绑定数据和元数据的能力。有三个主要功能，从Biobase包，提供对实验级元数据、主题级元数据和数据本身的访问。

来访问实验级元数据experimentData ()函数用于返回aMIAME(Minimum Information About a Microarray Experiment)对象。

实验数据(zeller.eset)
##实验数据##实验人员姓名:Zeller G, Tap J, Voigt AY, Sunagawa S, Kultima JR, Costea PI, Amiot A, Böhm J, Brunetti F, Habermann N, Hercog R, Koch M, Luciani A, Mende DR, Schneider MA, Schrotz-King P, Tournigand C, Tran Van Nhieu J, Yamada T, Zimmermann J, Benes V, Kloor M, Ulrich CM, von Knebel Doeberitz M, Sobhani I, Bork P ##实验室:德国海德堡欧洲分子生物学实验室结构与计算生物学单元。标题:粪便微生物群在结肠直肠癌早期检测中的潜力。## URL: ## PMIDs: 25432777 ## ##摘要:一个175字的摘要可用。使用“抽象”方法。##备注:##测序平台:## IlluminaHiSeq

来访问主题级元数据pData ()函数用于返回adata.frame包含学科水平的研究变量。

head(pData(zeller.eset))
# # # # subjectID body_site antibiotics_current_use CCIS27304052ST-3-0 fr - 001份粪便< NA > # # CCIS13047523ST-4-0 fr - 003份粪便< NA > # # CCIS15794887ST-4-0 fr - 024份粪便< NA > # # CCIS94603952ST-4-0 fr - 025份粪便< NA > # # CCIS74726977ST-3-0 fr - 026份粪便< NA > # # CCIS02856720ST-4-0 fr - 027份粪便< NA > # # # # study_condition疾病disease_subtype年龄CCIS27304052ST-3-0控制< NA > < NA > 52 # # CCIS13047523ST-4-0控制< NA > < NA > 70 # # CCIS15794887ST-4-0控制< NA > < NA > 37 # # CCIS94603952ST-4-0控制< NA >< NA > 80 # # CCIS74726977ST-3-0腺瘤腺瘤smalladenoma 66 # # CCIS02856720ST-4-0控制< NA > < NA > 74 # # age_category性别BMI国家non_westernized # # CCIS27304052ST-3-0成年女性20联邦铁路局没有# # CCIS13047523ST-4-0高级男22 # # CCIS15794887ST-4-0成年女性18联邦铁路局没有# # CCIS94603952ST-4-0高级女21联邦铁路局没有# # CCIS74726977ST-3-0高级男24 # # CCIS02856720ST-4-0高级男性27联邦铁路局没有# # DNA_extraction_kit number_reads number_bases # # CCIS27304052ST-3-0 Gnome48552680 NA # # CCIS13047523ST-4-0 Gnome 48552680 NA # # CCIS15794887ST-4-0 Gnome 45764894 NA # # CCIS94603952ST-4-0 Gnome 45764894 NA # # CCIS74726977ST-3-0 Gnome 140945980 NA # # CCIS02856720ST-4-0 Gnome 140945980 NA # # minimum_read_length median_read_length与fobt tnm 45 # # CCIS27304052ST-3-0 NA < NA >没有< NA > # # CCIS13047523ST-4-0 45 NA < NA >没有< NA > # # CCIS15794887ST-4-0 45 NA < NA >没有< NA > # # CCIS94603952ST-4-0 1 NA < NA >没有< NA > # # CCIS74726977ST-3-0 45 NA < NA > < NA > # # CCIS02856720ST-4-0 1 NA yes  ## NCBI_accession ## CCIS27304052ST-3-0 ERR480588;ERR480587;ERR479092;ERR479091 ## CCIS13047523ST-4-0 ERR480532;ERR480529;ERR480531;ERR480530;ERR479036;ERR479034;ERR479035;ERR479033 ## CCIS15794887ST-4-0 ERR480538;ERR480537;ERR479042;ERR479041 ## CCIS94603952ST-4-0 ERR480895;ERR480894;ERR480893;ERR480896;ERR479400;ERR479399;ERR479398;ERR479397 ## CCIS74726977ST-3-0 ERR480794;ERR479298 ## CCIS02856720ST-4-0ERR480469; ERR480468; ERR480467 ERR480466; ERR478973 ERR478972; ERR478971; ERR478970

要访问数据本身(在本例中为相对丰度)exprs ()函数返回一个变量按样本(行按列)的数值矩阵。请注意，在分类学的所有层次上都存在“合成”分支，从王国开始，例如k__bacteria:

Exprs (zeller.eset)[1:6, 1:5] #前6行和5列
## CCIS27304052ST-3-0 CCIS13047523ST-4-0 ## k__Bacteria 95.80672 98.39758 ## k__Bacteria|p__Bacteroidetes 21.37090 24.69795 ## k__Archaea|p__Euryarchaeota 3.11078 0.00847 ## CCIS15794887ST-4-0 ## k__Bacteria 99.99562 100.00000 ## k__Archaea 0.00438 0.00000 ## k__Viruses 0.00000 0.00000 ## k__Bacteria|p__Firmicutes 55.31008 71.41404 ## k__Bacteria|p__Bacteroidetes 23.7370616.01887 ## k__Archaea|p__Euryarchaeota 0.00438 0.00000 ## CCIS74726977ST-3-0 ## k__Bacteria 90.08853 ## k__Archaea 1.98404 ## k__Viruses 7.92743 ## k__Bacteria|p__Firmicutes 58.24313 ## k__Bacteria|p__Bacteroidetes 27.91643 ## k__Archaea|p__Euryarchaeota 1.98404

Bioconductor提供了ExpressionSet类的进一步文档，并发表了一篇优秀的介绍．

创建phyloseq和metagenomeSeq对象

curatedMetagenomicData提供方便的转换函数ExpressionSet对象phyloseq: phyloseq和metagenomeSeq: MRExperiment用于下游分析的类对象。例如，要转换zeller.eset对象指向这两个类:

ExpressionSet2phyloseq (zeller.eset)
## phyloseq-class实验级对象## otu_table() OTU表:[10503个分类单元，199个样本]## sample_data()样本数据:[199个样本，21个样本变量]## tax_table()分类表:[10503个分类单元，8个分类等级]
ExpressionSet2MRexperiment (zeller.eset)
## mexperiment (storageMode: environment) ## assayData: 10503个特征，199个样本##元素名称:计数## protocolData: none ## phenoData ## sampleNames: CCIS27304052ST-3-0 CCIS13047523ST-4-0…## MMPU72854103ST (199 total) ## varLabels: subjectID body_site…NCBI_accession (21 total) ## varMetadata: labelDescription ## featureData ## featureNames: k__Bacteria k__Archaea…t__PRJNA15006 (10503 ## total) ## fvarLabels: Kingdom Phylum…Strain(共8个)## fvarMetadata: labelDescription ## experimentData:使用'experimentData(object)' ##注释:

这也适用于合并eset从上面的对象:

ExpressionSet2phyloseq (eset)
## phyloseq-class实验级对象## otu_table() OTU表:[1914个分类群和172个样本]## sample_data()样本数据:[172个样本，19个样本变量]## tax_table()分类表:[1914个分类群，8个分类等级]
ExpressionSet2MRexperiment (eset)
## mexperiment (storageMode: environment) ## assayData: 1914个特征，172个样本##元素名称:计数## protocolData:无## phenoData ## sampleNames: BritoIL_2016.metaphlan_bugs_list.oralcavity:M1.1。## britoil_2016 .metaphlan_bugs_list.oral - cavity:M1.10。山……# # Castro-NallarE_2015.metaphlan_bugs_list。或alcavity:ES_080 (172 ## total) ## varLabels: body_site antibiotics_current_use ... studyID (19 ## total) ## varMetadata: labelDescription ## featureData ## featureNames: k__Bacteria k__Viruses ... ## t__Jonquetella_anthropi_unclassified (1914 total) ## fvarLabels: Kingdom Phylum ... Strain (8 total) ## fvarMetadata: labelDescription ## experimentData: use 'experimentData(object)' ## Annotation:

的例子phyloseq使用curatedMetagenomicData数据集的分析提供在curatedMetagenomicData装饰图案．

系统发育树和UniFrac距离

集phylogenetictree = TRUE在…中包含一个系统发育树phyloseq对象(仅限生物开发):

泽勒。<- ExpressionSet2phyloseq(zeller.eset, phylogenetictree = TRUE)
wt = UniFrac(zeller);tree, weighted=TRUE, normalized=FALSE, parallel=FALSE, fast=TRUE) plot(hclust(wt)， main=“加权UniFrac距离”)

大肠杆菌流行率的探索性分析

这里有一个直接的，探索性的分析大肠杆菌在zeller数据集中使用ExpressionSet对象。类提供的子集方法将在稍后提供更优雅的解决方案phyloseq包，但对于用户熟悉grep ()和ExpressionSet对象，这样的手动方法就足够了。

首先,这大肠杆菌有相关的分类群吗?

grep("coli"， rownames(zeller.eset)， value=TRUE)
##[1]“k_细菌|p_变形菌门|c_伽马变形菌门|o_肠杆菌门|f_肠杆菌科|g_埃希菌门|s_埃希菌大肠杆菌”##[2]”k_细菌|p_变形菌门|c_伽马变形菌门|o_肠杆菌门|f_肠杆菌科|g_埃希菌门|s_埃希菌门大肠杆菌|t_埃希菌门大肠杆菌未分类“##[3]”k_细菌|p_厚壁菌门|c_梭菌门|o_梭菌门|f_瘤胃球菌科|g_无氧主干|s_无氧主干大肠杆菌”## [4]"k_细菌|p_厚壁菌门|c_梭菌门|c_梭菌门|f_瘤胃球菌科|g_厌氧菌门|s_厌氧菌门|c_梭菌门|t_ gcf_000154565 " ## [5] "k_细菌|p_厚壁菌门|c_梭菌门|f_胃链球菌科|g_胃链球菌科|p_厚壁菌门|c_梭菌门|o_梭菌门|f_胃链球菌科|g_胃链球菌科|s_梭菌门糖colicum|t_ gcf_000373865 " ## [7]“k_细菌|p_变形菌门|c_ epsilon变形菌门|o_弯曲菌门|f_弯曲菌科|g_弯曲菌门|s_弯曲菌科|o_弯曲菌门|c_弯曲菌科|g_弯曲菌门|s_弯曲菌门大肠|t_弯曲菌门大肠未分类”##[9]“k_细菌|p_放线菌门|c_放线菌门|o_放线菌门|f_分枝杆菌科|g_ amycolicoccus”## [10]“k_细菌|p_厚壁菌门|c_梭菌门|o_梭菌门|f_胃球菌科|g_合滋养肉杆菌门|s_合滋养肉毒杆菌门|s_合滋养肉毒杆菌门|s_放线菌门|f_分枝杆菌科|g_ amycolicoccoccus |s_ amycolicoccus_subflavus”##[12]”k_细菌|p_厚壁菌门|c_梭菌门|o_梭菌门|f_胃球菌科|g_合滋养肉毒杆菌门|s_合滋养肉毒杆菌门糖结肠|t_ gcf_000190635”## [13]"k_细菌|p_放线菌门|c_放线菌门|o_放线菌门|f_分枝杆菌科|g_ amycolicoccus |s_ amycolicoccus_subflavus |t_ gcf_000214175 " ## [14] "k_细菌|p_变形菌门|c_ γ变形菌门|o_肠杆菌门|f_肠杆菌科|g_耶尔森氏菌|s_耶尔森氏菌小肠结肠炎" ## [15]"k_细菌|p_变形菌门|c_ γ变形菌门|o_弧菌门|f_弧菌科|g_弧菌|s_海带弧菌" ## [16]“k_细菌|p_螺旋体|c_螺旋体|c_螺旋体|f_短螺旋体|g_短螺旋体|s_短螺旋体|f_短螺旋体|g_短螺旋体|s_短螺旋体|f_肠杆菌科|g_耶尔森菌|s_耶尔森菌|c_耶尔森菌|o_弧菌|f_弧菌科|g_弧菌|s_海带弧菌|”## [19]"k_细菌|p_螺旋体纲|c_螺旋体纲|o_螺旋体纲|f_短螺旋体纲|g_短螺旋体纲|s_短螺旋体纲|s_短螺旋体纲|s_短螺旋体纲|c_厚壁菌纲|c_梭菌纲|g_梭菌纲|s_梭菌纲colicanis" ## [21] "k_细菌|p_厚壁菌纲|c_梭菌纲|o_梭菌纲|f_梭菌科|g_梭菌纲|s_梭菌纲colicanis|t__GCF_000371465"

创建一个向量大肠杆菌相对丰度。这grep调用的最后带" $ "选择唯一以" s__Escherichia_coli "结尾的行:

x = exps (zeller.eset) [grep("s__Escherichia_coli$"， rownames(zeller.eset))，] summary(x)
* * * * * * * * * * * * * * * *## 0.00000 0.01985 0.15273 1.64121 1.37392 37.54376

这可以绘制为直方图:

hist(x, xlab =“相对丰度”，main=“大肠杆菌患病率”，break =“FD”)

估计绝对原始计数数据

的列相乘，可以从相对计数数据中估计绝对原始计数数据ExpressionSet数据中每个样本的读取次数，如pData列“number_reads”。出于演示目的，你可以(但不是必须!)通过除以100再乘以读取次数来手动完成:

泽勒。计数= sweep(exprs( zeller.eset ), 2, zeller.eset$number_reads / 100, "*") zeller.counts = round(zeller.counts) zeller.counts[1:6, 1:5]
## CCIS27304052ST-3-0 CCIS13047523ST-4-0 ## k__Bacteria 46516730 61883930 ## k__Archaea 1510367 537791 45796391 ## k__Bacteria|p__Bacteroidetes 10376145 15532965 ## k__Archaea|p__Euryarchaeota 1510367 5327 ## CCIS15794887ST-4-0 ## k__Bacteria 45762889 69804753 ## k__Archaea 2005 0 ## k__Viruses 00 ## k__Bacteria|p__Firmicutes 25312599 49850394 ## k__Bacteria|p__Bacteroidetes 10863240 11181933 ##k__Archaea|p__Euryarchaeota 2005 0 ## CCIS74726977ST-3-0 ## k__Bacteria 126976161 ## k__Archaea 2796425 ## k__Viruses 11173394 ## k__Bacteria|p__Firmicutes 82091350 ## k__Bacteria|p__Bacteroidetes 39347086 ## k__Archaea|p__Euryarchaeota 2796425

或者直接设置计数论点curatedMetagenomicData ()来真正的：

zeller.eset2 = curatedMetagenomicData("ZellerG_2014.metaphlan_bugs_list. exe ")“大便”，counts = TRUE, dryrun = FALSE)[b[1]]
正在开发ZellerG_2014.metaphlan_bugs_list.stool
## snapshotDate(): 2017-06-09
请参阅curatedMetagenomicData和browseVignettes('curatedMetagenomicData')获取文档
##从缓存加载'/home/ubuntu//ExperimentHub / 535”
all.equal (exprs (zeller.eset2)、zeller.counts)
## [1] true

分类法感知分析使用phyloseq

对于MetaPhlAn2 bug数据集(但不包括其他数据类型)，通过转换为a，您可以获得许多系统发育感知的生态分析和绘图phyloseq类对象。curatedMetagenomicData提供了ExpressionSet2phyloseq ()函数使这变得简单:

泽勒suppressPackageStartupMessages(库(phyloseq))。seq = ExpressionSet2phyloseq(zeller.eset)

请注意下列对ExpressionSet2phyloseq有用的参数:

使用Metavizr进行可视化

现在我们可以使用Metavizr检查zeller数据集。这个可视化工具允许对分类层次和统计引导的可视化分析进行交互式探索。具体来说，metavizr通过一种名为FacetZoom的导航机制来实现功能选择。我们首先为zeller数据集创建一个mexperimentExpressionSet2MRexperiment．为了在本次研讨会中进行演示，我们将对60岁以下的样本进行分析，并对对照或CRC研究条件进行分析。

zeller_MR_expr = expressionset2mexperiment (zellerer .eset) zeller_MR_expr <- zeller_MR_expr[， which(pData(zeller_MR_expr)$age < 60)] zeller_MR_expr <- zeller_MR_expr[，-which(pData(zeller_MR_expr)$study_condition == "adenoma")] public_ipv4 <- try(httr::content(httr::GET("http://169.254.169.254/latest/meta-data/public-ipv4"))， silent = TRUE) if(grepl("Timeout is reached")，simpleMessage(public_ipv4)) {app <- startMetaviz(host="http://metaviz-dev.cbcb.umd.edu")} else{app <- startMetaviz(host="http://metaviz-dev.cbcb.umd.edu"， ws_host = public_ipv4, daemonized = FALSE, try_ports = TRUE)}

打开浏览器后，我们可以从zeller_MR_expr对象中可视化宏基因组特征。我们使用情节方法执行此操作。的情节函数可以接受在epivirData中注册的任何类作为输入。

生物导体AMI注:在rstudio-server版本< 1.1.273的情况下，用户需要在每次用户在UI上输入后运行app$service()命令。Rstudio-server有一个干扰R事件循环的bug。Metavizr应该可以在rstudio桌面上正常运行。

facetZoom <- app$plot(zeller_MR_expr, type = "InnerNodeCounts"， datasource_name = "zeller"， feature_order = colnames(fData(zeller_MR_expr))) app$service()

从元基因组图像中:您现在应该看到一个FacetZoom来探索来自mexperiment对象的元基因组特征的层次结构。为了导航复杂的、层次结构的特征空间，我们使用了FacetZoom可视化设计，并将其扩展到宏基因组数据。由于屏幕大小和渲染大树的性能的限制，这种技术是导航树的有效可视化。FacetZoom可以帮助放大和缩小树，并遍历子树。树中的每个节点都有一个与之相关的状态，并且节点有三种可能的状态1)扩展-在分析期间使用所有子树节点2)崩溃-将子树下的所有节点聚集到选定的节点3)删除-丢弃子树中用于分析的所有节点。节点的状态也会传播到它的所有子节点，并且可以通过节点的不透明度来识别。左侧提供行级控件—用于设置层次结构中选定深度/分类类的所有节点的状态。FacetZoom的右侧显示了当前子树的位置。用户可以在节点上设置状态来定义特征空间上的切割。切割定义了计数数据如何在链接的图和图表中可视化。 In addition to defining the cut, we also have a navigation bar to limit the range of features when querying count data. Navigation controls move the bar left/right and extend over the entire range of the current tree in the icicle. Each of these actions queries the underlying data structure and automatically propagates the changes to other visualizations in the workspace. When navigating outside the scope of the navigation bar, chevrons (left/right) appear on the navigation bar to help identify the current position. Hovering over a node highlights the entire lineage in the FacetZoom along with all other linked charts and plots.

现在我们可以添加一个热图，显示来自mexperiment的计数数据。后应用程序服务()美元调用，当R命令提示符上显示“Serving Epiviz, escape to continue interactive session…”时，我们可以点击节点标签的中心图标P在FacetZoom的左手边这将更新热图，从显示类级别的计数到显示门级别的计数。消息“zeller_1 datasourceGroup缓存已清除”应该出现在R会话中。我们也可以在UI中进行选择并运行应用程序服务()美元处理请求。

<- app$chart_mgr$revisualize(chart_type = "HeatmapPlot"， chart = facetZoom) app$service()

通过在FacetZoom中设置节点状态，我们可以在分类法的任何级别上删除血统并显示计数。一旦我们完成了对数据的研究，我们就可以通过调用图表管理器来删除这些图表。

应用chart_mgr rm_all_charts美元()

使用metagenomeSeq计算差分丰度，并在metavez中使用特征选择

现在我们将使用metagenomeSeq来计算zeller数据集在不同分类级别上的差异丰度。对于metagenomeSeq数据模型，我们将只保留叶计数，然后使用metagenomeSeq函数aggregateByTaxonomy在层次结构的不同级别上计算计数。

#删除除最低级别zeller_MR_expr之外的所有级别的分类法计数和特征zeller_MR_expr <- zeller_MR_expr[-unname(which(is.na(fData(zeller_MR_expr)$Strain))))，] #删除mexperimental计数数据rownames(zeller_MR_expr) <- sapply(strsplit(rownames(zeller_MR_expr)， "__")， function(i){unname(i)[2]})中每个叶级特征的t__

现在，我们可以使用fitFeatureModel函数来估计分类中最低级别的特征(在本例中为Strain)在样本之间的log fold变化控制和儿童权利公约研究条件。我们首先使用累积和缩放方法对计数数据进行归一化。

#使用累积和缩放方法规范化计数。zeller_MR_expr <- filterData(zeller_MR_expr, present = 5) zeller_MR_expr <- cumNorm(zeller_MR_expr, p = 0.75) #设置模型为study_condition，并使用fitFeatureModel测试zeller_sample_data <- pData(zeller_MR_expr) mod <- model。矩阵(~1+study_condition, data = zeller_sample_data) results_zeller <- fitFeatureModel(zeller_MR_expr, mod) logFC_zeller <- MRcoefs(results_zeller, number = nrow(results_zeller))

metagenomeSeq提供了多种测试方法。我们还可以使用fitZig方法来测试study_condition对性别的控制。fitZig将后验概率添加到mexperiment中，因此我们创建了一个单独的mexperiment来传递给该函数。

#创建一个单独的MRexperiment传给fitZig zig_zeller_MR_expr = ExpressionSet2MRexperiment (zeller.eset) zig_zeller_MR_expr < - zig_zeller_MR_expr [, pData (zig_zeller_MR_expr)年龄< 60美元)]zig_zeller_MR_expr < - zig_zeller_MR_expr[,——(pData (zig_zeller_MR_expr) $ study_condition = =“腺瘤”)]zig_zeller_MR_expr < - zig_zeller_MR_expr [-unname ((is.na (fData (zig_zeller_MR_expr)应变美元))),]rownames (zig_zeller_MR_expr) < -酸式焦磷酸钠(strsplit (rownames (zig_zeller_MR_expr),“__”),function(i){unname(i)[2]}) zig_zeller_MR_expr <- filterData(zig_zeller_MR_expr, present = 5) zig_zeller_MR_expr <- cumNorm(zig_zeller_MR_expr, p = 0.75) zig_zeller_sample_data <- pData(zig_zeller_MR_expr) mod_zig <- model。矩阵(~1+study_condition+gender, data = zig_zeller_sample_data) results_zeller_zig <- fitZig(zig_zeller_MR_expr, mod_zig) coefs_zeller_zig <- MRtable(results_zeller_zig, number = nrow(results_zeller_zig))

metagenomeSeq包插图提供了关于实现的可用测试以及对分析有用的其他功能的进一步详细信息。

我们将使用fitFeatureModel的结果，因为它为控制组和crc组中样本之间的这个级别的特征提供了logFC估计。我们可以用表格的形式来检查它们。我们还可以使用这些结果来选择在Metaviz中不需要考虑的特性名称。在使用logFC阈值和调整后的p值截止值删除特征后，我们可以交互式地探索不同层次结构中剩余特征的结果。

features <- rownames(logFC_zeller) featuresToKeep_names <- features[which(logFC_zeller[which(abs(logFC_zeller$logFC) > 2)，]$adjPvalues <.1)] featuresToKeep <- rep(2, length(featuresToKeep_names)) names(featuresToKeep) <- featuresToKeep_names featuresToRemove_names <- features[!](features %in% featustokeep_names)] featuresToRemove <- rep(0, length(featuresToRemove_names)) names(featuresToRemove) <- featuresToRemove_names featureSelection = c(featustokeep, featuresToRemove)

我们可以添加新的FacetZoom、热图和堆叠图，然后选择感兴趣的特征。一旦添加了热图和堆叠图，我们就可以导航到分类法的最低级别，并查看从考虑中删除的特性。

facetZoom <- app$plot(zeller_MR_expr, type =" LeafCounts"， datasource_name="zeller_leaf_level"， feature_order = colnames(fData(zeller_MR_expr))) app$service()热图<- app$chart_mgr$revisualize(chart_type =" HeatmapPlot"， chart = facetZoom) app$service() stackedPlot <- app$chart_mgr$revisualize(chart_type ="StackedLinePlot"， chart = facetZoom) app$service()

现在我们可以调用propagateHierarchyChanges具有趣味性的特点。应删除不显示差异丰度的应变特征。

app$get_ms_object(facetZoom)$propagateHierarchyChanges(featuresselection, request_with_labels = TRUE)

我们还可以使用聚合到分类法另一层的计数执行相同的fitFeatureModel计算，并将这些更改传播到现有的可视化中。在本例中，我们将计数聚合到分类法的Class级别，然后只保留在两组之间显示logFC大于1的特性。

aggregation_level = "Class" agg_zeller_MR_expr <- aggregateByTaxonomy(zeller_MR_expr, lvl=aggregation_level) agg_zeller_MR_expr <- cumNorm(agg_zeller_MR_expr, p = 0.75) agg_mod <- model。matrix(~1+study_condition, data = pData(agg_zeller_MR_expr)) agg_results_zeller <- fitFeatureModel(agg_zeller_MR_expr, agg_mod) agg_logFC_zeller <- MRcoefs(agg_results_zeller, number = nrow(agg_results_zeller)) agg_features <- rownames(agg_logFC_zeller) agg_featuresToKeep_names <- agg_features[which(agg_logFC_zeller[which(abs(agg_logFC_zeller$logFC) > 1)，]$adjPvalues <.1)] agg_featuresToKeep <- rep(2，length(agg_featuresToKeep_names)) names(agg_featuresToKeep) <- agg_featuresToKeep_names agg_featuresToRemove_names <- agg_features[!](agg_features %in% agg_featuresToKeep_names)] agg_featuresToRemove <- rep(0, length(agg_featuresToRemove_names)) names(agg_featuresToRemove) <- agg_featuresToRemove_names agg_selectionUpdate <- c(agg_featuresToKeep, agg_featuresToRemove) app$get_ms_object(facetZoom)$propagateHierarchyChanges(agg_selectionUpdate, request_with_labels = TRUE)