README.md
In mengqingren/VTwins:
VTwins
The goal of VTwins is to perform phenotype-enriched features including species or pathways for metagenomic sequencing or 16S sequencing.
if you want to use standalone version in Linux, you can check VTwins.Linux from mengqingren/VTwins.Linux.
You can also find the R script of this paper mengqingren/VTwins.Linux.
if you want to use R package, you can keep reading.
Installation
Now VTwins is not on cran, You can install the development version of
VTwins from GitHub with:
``` r
install.packages("remotes")
remotes::install_github("mengqingren/VTwins")
### Dependence
- R version 4.0.5
- tidyverse_1.3.1
### Basic Usage
``` r
library(VTwins)
pair_find(data=YourRelativeAbundanceDataframe, # must be a data frame with columns representing features and rows representing samples.
phenodata=YourPhenotypeDataframe, # must be a data frame with two columns, and colnames are `id` (column 1) and `grp` (column 2). column id represent the sample id, column grp consist of `grp1` and `grp2`, representing the ctrl and disease, repsectively.
k="euclidean", # distance calculating method. it must bu consist with the method in `dist` function.
SavePath = NULL, # output directory. Default: ./
Cut_pair=25,
method_choose=c("Wilcox","Permutation"),
ShuffleWstat = NULL, # output W stat
BoundarySample = NULL, # output boundary samples with distance
BoundaryPair=NULL, # output final pairs
ShuffleTime=10000, # shuffle time
DownPercent = 0.2, # lower percentage of shuffle pairs
Uppercent=0.8,
PvalueCutoff = 0.05) # p value cutoff of `Incre.aveRank.P` and `Decre.aveRank.P`
Parameter Description
- data: must be a data frame with columns representing features and rows representing samples.
- phenodata: must be a data frame with two columns, and colnames are
id (column 1) and grp (column 2). column id represent the sample id, column grp consist of grp1 and grp2, representing the ctrl and disease, repsectively.
- k: distance calculating method. it must bu consist with the method in
dist function.
- SavePath: filename of output directory. Default: ./
- Cut_pair: the cutoff of redundant pairs to perform permutation or wilcox rank paired test. Defult 25
- method_choose: necessary parameter of greated than 25 pair, choose from :"Wilcox","Permutation",
- ShuffleWstat: filename of shuffle W stats
- BoundarySample: filename of output boundary samples with distance
- BoundaryPair: filename of output final pairs
- ShuffleTime: shuffle time
- DownPercent: lower percentage of shuffle pairs
- Uppercent: higher percentage of shuffle pairs
- PvalueCutoff: p value cutoff of
Incre.aveRank.P and Decre.aveRank.P
Note: if the sample pairs are less than 10, it will return nothing.
Output
- Decre.aveRank.P : P value based on averange rank. Generally, to confirm the significant features, we usually use variable
Decre.aveRank.P to evaluate the disease-enriched features and Incre.aveRank.P to evaluate the control enriched features.
- Incre.aveRank.P : P value based on averange rank. Generally, to confirm the significant features, we usually use variable
Decre.aveRank.P to evaluate the disease-enriched features and Incre.aveRank.P to evaluate the control enriched features.
- Decre.minRank.* : P value based on min rank.
- Incre.minRank.P : P value based on min rank.
- Decre.maxRank.P : P value based on max rank.
- Incre.maxRank.P : P value based on max rank.
- Species : Feature
- Enriched : phenotype-enriched groups with provided p value cutoff
- Decre.maxRank.P.FDR : P.adjust of Decre.maxRank.P
- Decre.minRank.P.FDR : P.adjust of Decre.minRank.P
- Decre.aveRank.P.FDR : P.adjust of Decre.aveRank.P
- Incre.aveRank.P.FDR : P.adjust of Incre.aveRank.P
- Decreasing.Rank.Max : shuffle W stats rank based on max method by decreasing
- Increasing.Rank.Max : shuffle W stats rank based on max method by increasing
- Decreasing.Rank.Min : shuffle W stats rank based on min method by decreasing
- Increasing.Rank.Min : shuffle W stats rank based on min method by increasing
- Decreasing.Rank.Average : shuffle W stats rank based on average method by decreasing
- Increasing.Rank.Average : shuffle W stats rank based on average method by increasing
- Ctlmean: mean value of relativa abundance for paired control samples
- Dismean: mean value of relativa abundance for paired disease samples
Important : Must keep the Same Order of samples in relative abundance dataframe and phenotype dataframe. And we also keep the samples clustered and placed according the grp1 and grp2 of variable grp in phenotype data.
You can refer to the following example data format.
Usage -> Mock Data
if you want to skip the download step, you can find the mock and real datasets in mengqingren/VTwins.Linux
library(tidyverse)
library(VTwins)
library(vegan)
# relative abundance dataframe
set.seed(12345)
dataset <- data.frame(matrix(runif(1200, min = 1e-5, max = 1),nrow = 120,ncol = 10))
colnames(dataset) <- paste("Feature",1:10,sep = '')
rownames(dataset) <- paste("Sample",1:120,sep = '')
dataset.normalized <- decostand(dataset,method = "total",1) #normalization for fetures like species's relative abundance
#write.table(dataset.normalized,file = "test.data.txt",sep = '\t',quote = F)
# phenotype dataframe
phe_data <- data.frame(id = paste("Sample",1:120,sep = ''),grp=rep(c("grp1","grp2"),c(60,60)))
### dataset.normalized <- read.table("test.data.txt",header = T,row.names = 1,sep = '\t')
### phe_data <- read.table("test.phenodata.txt",header = T,sep = "\t")
# Run
res <- pair_find(data=dataset.normalized,
phenodata=phe_data,
k="euclidean",
Cut_pair=25,
method_choose="Permutation",
SavePath = "./",
ShuffleWstat = "ShuffleWstat",
BoundarySample = "BoundarySample",
BoundaryPair="BoundaryPair",
ShuffleTime=10000,
DownPercent = 0.2,
Uppercent=0.8,PvalueCutoff=0.01)
res$log2FC = log2(as.numeric(res$Dismean)/as.numeric(res$Ctlmean))
write.csv(res,"Results.csv",row.names=F)
Usage -> Real Data
library(curatedMetagenomicData) #bioconductor
library(phyloseq) #bioconductor
library(tidyverse)
# download dataset -> object
ZellerGData <- curatedMetagenomicData("ZellerG_2014.metaphlan_bugs_list*",
dryrun = FALSE,
counts = TRUE,
bugs.as.phyloseq = TRUE)
psAB <- subset_samples(ZellerGData$ZellerG_2014.metaphlan_bugs_list.stool, study_condition != "adenoma")
psAB <- prune_samples(sample_sums(psAB) >= 10^3, psAB)
CRC_WMS <- filter_taxa(psAB,function(x) sum(x>0)>0,1)
# phenodata frame
grp1_name = "control"
grp2_name = "CRC"
variable_name = "study_condition"
# transitional phenotype dataframe
phe_data <- CRC_WMS@sam_data %>% data.frame() %>% dplyr::select("study_condition") %>% rownames_to_column() %>% dplyr::rename(id=1,grp=2) %>%
mutate(grp = dplyr::case_when(grp == "control" ~ "grp1",grp == "CRC" ~ "grp2")) %>% arrange(grp)
# transitional relative abundance dataframe
dataset <- t(CRC_WMS@otu_table@.Data)%>% data.frame() %>% rownames_to_column() %>% merge(phe_data,.,by.x="id",by.y="rowname") %>%
arrange(grp)
# Final phenotype dataframe
phe_data2 <- dataset %>% dplyr::select(id,grp)
# Final relative abundance dataframe
dataset2 <- dataset %>% dplyr::select(-grp) %>% remove_rownames() %>% column_to_rownames("id")
### dataset2 <- readRDS("RealData.dataset2.RDS")
### phe_data2 <- readRDS("RealData.phe_data2.RDS")
# Run
res <- pair_find(data=dataset2,
phenodata=phe_data2,
k="euclidean",
Cut_pair=25,
method_choose="Permutation",
SavePath = "./",
ShuffleWstat = "RealData.ShuffleWstat",
BoundarySample = "RealData.BoundarySample",
BoundaryPair="RealData.BoundaryPair",
ShuffleTime=10000,
DownPercent = 0.2,
Uppercent=0.8,PvalueCutoff=0.01)
res$log2FC = log2(as.numeric(res$Dismean)/as.numeric(res$Ctlmean))
write.csv(res,"RealData.Results.csv",row.names=F)
mengqingren/VTwins documentation built on July 1, 2023, 9:16 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
VTwins
The goal of VTwins is to perform phenotype-enriched features including species or pathways for metagenomic sequencing or 16S sequencing.
if you want to use standalone version in Linux, you can check VTwins.Linux from mengqingren/VTwins.Linux.
You can also find the R script of this paper mengqingren/VTwins.Linux.
if you want to use R package, you can keep reading.
Installation
Now VTwins is not on cran, You can install the development version of
VTwins from GitHub with:
``` r
install.packages("remotes")
remotes::install_github("mengqingren/VTwins")
### Dependence
- R version 4.0.5
- tidyverse_1.3.1
### Basic Usage
``` r
library(VTwins)
pair_find(data=YourRelativeAbundanceDataframe, # must be a data frame with columns representing features and rows representing samples.
phenodata=YourPhenotypeDataframe, # must be a data frame with two columns, and colnames are `id` (column 1) and `grp` (column 2). column id represent the sample id, column grp consist of `grp1` and `grp2`, representing the ctrl and disease, repsectively.
k="euclidean", # distance calculating method. it must bu consist with the method in `dist` function.
SavePath = NULL, # output directory. Default: ./
Cut_pair=25,
method_choose=c("Wilcox","Permutation"),
ShuffleWstat = NULL, # output W stat
BoundarySample = NULL, # output boundary samples with distance
BoundaryPair=NULL, # output final pairs
ShuffleTime=10000, # shuffle time
DownPercent = 0.2, # lower percentage of shuffle pairs
Uppercent=0.8,
PvalueCutoff = 0.05) # p value cutoff of `Incre.aveRank.P` and `Decre.aveRank.P`
Parameter Description
- data: must be a data frame with columns representing features and rows representing samples.
- phenodata: must be a data frame with two columns, and colnames are
id(column 1) andgrp(column 2). column id represent the sample id, column grp consist ofgrp1andgrp2, representing the ctrl and disease, repsectively. - k: distance calculating method. it must bu consist with the method in
distfunction. - SavePath: filename of output directory. Default: ./
- Cut_pair: the cutoff of redundant pairs to perform permutation or wilcox rank paired test. Defult 25
- method_choose: necessary parameter of greated than 25 pair, choose from :"Wilcox","Permutation",
- ShuffleWstat: filename of shuffle W stats
- BoundarySample: filename of output boundary samples with distance
- BoundaryPair: filename of output final pairs
- ShuffleTime: shuffle time
- DownPercent: lower percentage of shuffle pairs
- Uppercent: higher percentage of shuffle pairs
- PvalueCutoff: p value cutoff of
Incre.aveRank.PandDecre.aveRank.P
Note: if the sample pairs are less than 10, it will return nothing.
Output
- Decre.aveRank.P : P value based on averange rank. Generally, to confirm the significant features, we usually use variable
Decre.aveRank.Pto evaluate the disease-enriched features andIncre.aveRank.Pto evaluate the control enriched features. - Incre.aveRank.P : P value based on averange rank. Generally, to confirm the significant features, we usually use variable
Decre.aveRank.Pto evaluate the disease-enriched features andIncre.aveRank.Pto evaluate the control enriched features. - Decre.minRank.* : P value based on min rank.
- Incre.minRank.P : P value based on min rank.
- Decre.maxRank.P : P value based on max rank.
- Incre.maxRank.P : P value based on max rank.
- Species : Feature
- Enriched : phenotype-enriched groups with provided p value cutoff
- Decre.maxRank.P.FDR : P.adjust of Decre.maxRank.P
- Decre.minRank.P.FDR : P.adjust of Decre.minRank.P
- Decre.aveRank.P.FDR : P.adjust of Decre.aveRank.P
- Incre.aveRank.P.FDR : P.adjust of Incre.aveRank.P
- Decreasing.Rank.Max : shuffle W stats rank based on max method by decreasing
- Increasing.Rank.Max : shuffle W stats rank based on max method by increasing
- Decreasing.Rank.Min : shuffle W stats rank based on min method by decreasing
- Increasing.Rank.Min : shuffle W stats rank based on min method by increasing
- Decreasing.Rank.Average : shuffle W stats rank based on average method by decreasing
- Increasing.Rank.Average : shuffle W stats rank based on average method by increasing
- Ctlmean: mean value of relativa abundance for paired control samples
- Dismean: mean value of relativa abundance for paired disease samples
Important : Must keep the Same Order of samples in relative abundance dataframe and phenotype dataframe. And we also keep the samples clustered and placed according the grp1 and grp2 of variable grp in phenotype data.
You can refer to the following example data format.
Usage -> Mock Data
if you want to skip the download step, you can find the mock and real datasets in mengqingren/VTwins.Linux
library(tidyverse)
library(VTwins)
library(vegan)
# relative abundance dataframe
set.seed(12345)
dataset <- data.frame(matrix(runif(1200, min = 1e-5, max = 1),nrow = 120,ncol = 10))
colnames(dataset) <- paste("Feature",1:10,sep = '')
rownames(dataset) <- paste("Sample",1:120,sep = '')
dataset.normalized <- decostand(dataset,method = "total",1) #normalization for fetures like species's relative abundance
#write.table(dataset.normalized,file = "test.data.txt",sep = '\t',quote = F)
# phenotype dataframe
phe_data <- data.frame(id = paste("Sample",1:120,sep = ''),grp=rep(c("grp1","grp2"),c(60,60)))
### dataset.normalized <- read.table("test.data.txt",header = T,row.names = 1,sep = '\t')
### phe_data <- read.table("test.phenodata.txt",header = T,sep = "\t")
# Run
res <- pair_find(data=dataset.normalized,
phenodata=phe_data,
k="euclidean",
Cut_pair=25,
method_choose="Permutation",
SavePath = "./",
ShuffleWstat = "ShuffleWstat",
BoundarySample = "BoundarySample",
BoundaryPair="BoundaryPair",
ShuffleTime=10000,
DownPercent = 0.2,
Uppercent=0.8,PvalueCutoff=0.01)
res$log2FC = log2(as.numeric(res$Dismean)/as.numeric(res$Ctlmean))
write.csv(res,"Results.csv",row.names=F)
Usage -> Real Data
library(curatedMetagenomicData) #bioconductor
library(phyloseq) #bioconductor
library(tidyverse)
# download dataset -> object
ZellerGData <- curatedMetagenomicData("ZellerG_2014.metaphlan_bugs_list*",
dryrun = FALSE,
counts = TRUE,
bugs.as.phyloseq = TRUE)
psAB <- subset_samples(ZellerGData$ZellerG_2014.metaphlan_bugs_list.stool, study_condition != "adenoma")
psAB <- prune_samples(sample_sums(psAB) >= 10^3, psAB)
CRC_WMS <- filter_taxa(psAB,function(x) sum(x>0)>0,1)
# phenodata frame
grp1_name = "control"
grp2_name = "CRC"
variable_name = "study_condition"
# transitional phenotype dataframe
phe_data <- CRC_WMS@sam_data %>% data.frame() %>% dplyr::select("study_condition") %>% rownames_to_column() %>% dplyr::rename(id=1,grp=2) %>%
mutate(grp = dplyr::case_when(grp == "control" ~ "grp1",grp == "CRC" ~ "grp2")) %>% arrange(grp)
# transitional relative abundance dataframe
dataset <- t(CRC_WMS@otu_table@.Data)%>% data.frame() %>% rownames_to_column() %>% merge(phe_data,.,by.x="id",by.y="rowname") %>%
arrange(grp)
# Final phenotype dataframe
phe_data2 <- dataset %>% dplyr::select(id,grp)
# Final relative abundance dataframe
dataset2 <- dataset %>% dplyr::select(-grp) %>% remove_rownames() %>% column_to_rownames("id")
### dataset2 <- readRDS("RealData.dataset2.RDS")
### phe_data2 <- readRDS("RealData.phe_data2.RDS")
# Run
res <- pair_find(data=dataset2,
phenodata=phe_data2,
k="euclidean",
Cut_pair=25,
method_choose="Permutation",
SavePath = "./",
ShuffleWstat = "RealData.ShuffleWstat",
BoundarySample = "RealData.BoundarySample",
BoundaryPair="RealData.BoundaryPair",
ShuffleTime=10000,
DownPercent = 0.2,
Uppercent=0.8,PvalueCutoff=0.01)
res$log2FC = log2(as.numeric(res$Dismean)/as.numeric(res$Ctlmean))
write.csv(res,"RealData.Results.csv",row.names=F)
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.