Phyloseq:表和元数据的相对丰富度不匹配

Question

我是一个相对较新的phyloseq，我努力获得一个相对丰富的siamcat R code for meta-analysis可接受的输入表。

# this works: from qza to phyloseq object
ps<-qza_to_phyloseq(
  features="all-table.qza",
  tree="rooted-tree.qza",
  taxonomy = "all-taxonomy.qza",
  metadata = "metafinal.tsv"
)

# import metadata
metadata <- read_tsv("metafinal.tsv")

# 30 overlap of the metadata-sample_id with ps, 115 only in metadata
gplots::venn(list(metadata=metadata$sample_id, features=sample_names(ps))

# works: from phyloseq object to relative abundance otu table
table(tax_table(ps)[, "Phylum"])
ps_rel_abund <- transform_sample_counts(ps, function(x){x / sum(x)})
ps_phylum_rel <- tax_glom(ps_rel_abund, "Phylum")
taxa_names(ps_phylum_rel) <- tax_table(ps_phylum_rel)[, "Phylum"]
rel_table <- as(otu_table(ps_phylum_rel), "matrix")

# column names and sample_id are 100% the same
colnames(rel_table)
metadata$sample_id

# 100% overlap:
gplots::venn(list(metadata=metadata$sample_id, featuretable=colnames(rel_table)))

# check that metadata and feature agree
stopifnot(all(colnames(rel_table) == metadata$sample_id))

在这里我得到了一个错误消息: all(colnames(rel_table) == metadata$sample_id)不是真的，并且下面的siamcat code根本不工作。

my metadata[1:5, 1:5]：sample_id absolute_filepath study experiment_acce…study_title

1 SRR8547628 $PWD/Chen_2020_da…Chen…c…的SRX5349649剖析2 SRR8547629 $PWD/Chen_2020_da…Chen…c…的SRX5349648剖析3 SRR8547630 $PWD/Chen_2020_da…Chen…c…的SRX5349647剖析4 SRR8547631 $PWD/Chen_2020_da…Chen…c…的SRX5349646剖析5 SRR8547632 $PWD/Chen_2020_da…Chen…c…的SRX5349645剖析

my rel-table[1:5, 1:5]：SRR5092146 SRR5092147 SRR5092148 SRR5092149芦荟0 0.0000000 0.00000000 0.000000000脊椎动物0 0.0000000 0.00000000 0.000000000复合丛枝菌0 0.0000000 0.00000000 0.000000000子囊菌0 0.0000000 0.00000000 0.000000000弯曲菌门0 0.2465222 0.01166882 0.004337051 SRR5092150鞭毛藻0.00000000脊椎动物0.00000000复合丛枝菌0.00000000子囊菌0.00000000弯曲菌门0.02106281

nrow(metadata)= 154 ncol(rel_table)= 154

为什么它不起作用呢？我尝试了几个星期，但我不能让代码正常运行……

感谢您的时间和帮助。

Answer 1

如果我正确理解了您的问题，您可能想知道，为什么元数据和要素表中的样本in之间有完美的重叠，但为什么

stopifnot(all(colnames(rel_table) == metadata$sample_id))

返回FALSE。

我想这是因为你的样品似乎顺序不同。元数据中的前五个示例是：

SRR8547628
SRR8547629
SRR8547630
SRR8547631
SRR8547632

您的特征表中的前五个示例是：

SRR5092146
SRR5092147
SRR5092148
SRR5092149
SRR5092150

试一试

stopifnot(all(colnames(rel_table) %in% metadata$sample_id))