本节主要是在PhyloseqTutorials上学习 (Pre)Processing Data

---

library(phyloseq)
#一般情况下让别人重复自己的分析做好加入软件的版本号
packagesVersion(phyloseq)
#加载包含在phyloseq包中的GlobalPatterns数据集。
data(GlobalPatterns)

Accessors

---

一个phyloseq对象的组成，比如OTU表，可以通过特殊的访问器函数或"Accessors"来访问。如果存在,于系统发育序列数据的特定信息,这些访问函数可用于用户和相关函数/包的直接交互.这些访问函数主要包括以下几个。

ntaxa()
#查看样品数
nsamples(GlobalPatterns)
sample_names(GlobalPatterns)[1:5]
#rank.names()和rank_names()有区别，rank.names()是去冗余的功能
rank_names(GlobalPatterns)
#获取sample_data中出现的示例变量
sample_variables(GlobalPatterns)
otu_table(GlobalPatterns)[1:5, 1:5]
tax_table(GlobalPatterns)[1:5, 1:4]
phy_tree(GlobalPatterns)
taxa_names(GlobalPatterns)[1:10]

myTaxa = names(sort(taxa_sums(GlobalPatterns), decreasing = TRUE)[1:10])
ex1 = prune_taxa(myTaxa, GlobalPatterns)
#注意plot_tree(ex1)和plot(phy_tree(ex1), show.node.label = TRUE)的区别
plot(phy_tree(ex1), show.node.label = TRUE)

Preprocessing

---

phyloseq包还包括筛选、设置和合并丰度数据的功能，在phyloseq中，Filtering采用了一种类似于genefilter包中的方法的模块化设计方式,这包括prune_taxa和prune_samples方法直接删除不需要的index,以及filterfun_sample和genefilter_sample函数，用于构建任意复杂的样例筛选条件,以及filter_taxa函数，用于按类别进行筛选,在下面的示例中，首先将GlobalPatterns数据转换为相对丰富度，创建新的GPr对象，然后对该对象进行过滤，使其只保留平均值大于10^-5的otu

#transform_sample_counts()该函数根据用户提供的函数转换一个分类单元丰度矩阵的样本计数
#function()为自定义函数
GPr  = transform_sample_counts(GlobalPatterns, function(x) x / sum(x) )
GPfr = filter_taxa(GPr, function(x) mean(x) > 1e-5, TRUE)

上述这就产生了一个高度细分的对象GPfr，它只包含原始的约19216个otu中的4624个。注意，在这两行中，我们分别提供了用于转换和过滤的自定义函数。
构建子集的方法prune_taxa和prune_samples适用于需要的otu或样本的完整子集直接可用的情况,另外，subset_taxa和subset_samples函数分别用于根据分类法表或Sample数据组件中包含的辅助数据进行子集设置。这些函数类似于core R中的subset其中初始数据参数后面是一个表示要保留的元素或行的任意逻辑表达式,此，根据有关辅助数据的条件表达式，可以将整个实验级数据对象作为子集。例如，下面的代码将首先分配给GP.chl是属于衣原体门的GlobalPatterns数据集的子集，然后删除总读取次数少于20次的样本

#主要的是区别prune_taxa()和subset_taxa()的差异
GP.chl = subset_taxa(GlobalPatterns, Phylum=="Chlamydiae")
GP.chl = prune_samples(sample_sums(GP.chl)>=20, GP.chl)

合并的方法包括merge_taxa和merge_samples分别用于合并特定的otu或样本。还有一个merge_phyloseq函数，用于两个或多个phyloseq对象(或一个phyloseq对象和一个或多个独立组件)的完全合并.例如，下面的代码合并了只包含衣原体的数据集中的前5个otu

GP.chl.merged = merge_taxa(GP.chl, taxa_names(GP.chl)[1:5])

在merge_taxa方法的基础上，phyloseq包包含了聚集函数，tip_glom和tax_glom,用于在一个实验中合并所有在系统发育或分类阈值之外相似的otu,下面的代码演示了如何在Family的分类级别上聚集“Bacteroidetes-only”数据集(称为gpsfb)，并创建结果的注释树

gpsfbg = tax_glom(gpsfb, "Family")
plot_tree(gpsfbg, color="SampleType", shape="Class", size="abundance")

对于通过任意R函数转换丰度值，phyloseq包含transform_sample_counts函数,它接受一个phyloseq对象和一个R函数作为参数，并返回一个phyloseq对象，其中的丰度值已经根据函数指定的转换进行了抽样转换,例如，下面的命令转换GP.chl丰度计数到百分数丰度

transform_sample_counts(GP.chl, function(OTU) OTU/sum(OTU) )

最后，以下是在创建主要的phyloseq手稿中的数字之前，应用于GlobalPatterns OTU计数的其余预处理步骤,移除至少20%的样本中不超过3次的分类单元,This protects against an OTU with small mean & trivially large C.V.

GP = filter_taxa(GlobalPatterns, function(x) sum(x > 3) > (0.2*length(x)), TRUE)

定义人与非人的分类变量，并将此新变量添加到样例数据中

sample_data(GP)$human = factor( get_variable(GP, "SampleType") %in% c("Feces", "Mock", "Skin", "Tongue") )

将丰度标准化到中位排序深度

#median为计算样本中位数
total = median(sample_sums(GP))
standf = function(x, t=total) round(t * (x / sum(x)))
gps = transform_sample_counts(GP, standf)

使用变异系数3.0的截止值筛选分类单元

#filter_taxa基于跨样本OTU丰度准则筛选分类单元
gpsf=filter_taxa(gps,function(x) sd(x)/mean(x)>3.0,TRUE)

提取拟杆菌的数据子集，用在一些绘图中

gpsfb = subset_taxa(gpsf, Phylum=="Bacteroidetes")

graphic summary

用图形来描述上述的这些数据

title = "plot_bar; Bacteroidetes-only"
plot_bar(gpsfb, "SampleType", "Abundance", title=title)
plot_bar(gpsfb, "SampleType", "Abundance", "Family", title=title)

其实本章节需要搞懂prune_taxa和filter_taxa以及subset_taxa的区别

• prune_taxa它的定义为从一个系统发育的对象中删除不需要的OTUs / taxa,例如从OTU表中过滤掉丰度小于10的OTUs/Taxa.prune_taxa(taxa,x)
• filter_taxa基于跨样本OTU丰度准则筛选分类单元,主要用于过滤来自phyloseq对象的otu，它接受一个phyloseq对象作为输入，并返回一个逻辑向量，该逻辑向量指示每个OTU是否通过了标准.filter_taxa(physeq, flist, prune=FALSE)
• subset_taxa从taxonomic expression表中提取物种记住这一点很重要，因为用户通常不知道这个子设置步骤还删除/省略了在表达式中某个地方有缺失值结果NA的OTUs，ex3 <- subset_taxa(GlobalPatterns, Phylum=="Bacteroidetes")