## ----echo=FALSE, results="hide", warning=FALSE,tidy=TRUE---------------------- suppressPackageStartupMessages({ library('LinkHD') }) ## ----results=' hide',eval=FALSE,tidy=TRUE------------------------------------- # #if you don't have devtools, please install it: # install.packages("devtools") # #if devtools is already installed, ommit previous step # library(devtools) # devtools::install_github(repo ="lauzingaretti/LinkHD" ,dependencies = TRUE) # library(LinkHD) ## ----results=' hide',eval=FALSE,tidy=TRUE------------------------------------- # ?? LinKData() ## ----message=FALSE,tidy=TRUE-------------------------------------------------- library(LinkHD) ## ----tidy=TRUE---------------------------------------------------------------- data("Taraoceans") data(Taraoceans) pro.phylo <- Taraoceans$taxonomy[ ,"Phylum"] TaraOc<-list(Taraoceans$phychem,as.data.frame(Taraoceans$pro.phylo), as.data.frame(Taraoceans$pro.NOGs)) TaraOc_1<-scale(TaraOc[[1]]) Normalization<-lapply(list(TaraOc[[2]],TaraOc[[3]]), function(x){DataProcessing(x,Method="Compositional")}) colnames(Normalization[[1]])=pro.phylo colnames(Normalization[[2]])=Taraoceans$GO TaraOc<-list(TaraOc_1,Normalization[[1]],Normalization[[2]]) names(TaraOc)<-c("phychem","pro_phylo","pro_NOGs") TaraOc<-lapply(TaraOc,as.data.frame) TaraOc<-lapply(TaraOc,as.data.frame) Tara_Out<-LinkData(TaraOc,Scale =FALSE,Center=FALSE,Distance = c("ScalarProduct","Euclidean","Euclidean")) ## ----tidy=TRUE---------------------------------------------------------------- CompromisePlot(Tara_Out)+theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), panel.background = element_blank(), axis.line = element_line(colour = "black")) ## ----tidy=TRUE---------------------------------------------------------------- CorrelationPlot(Tara_Out) + theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), panel.background = element_blank(), axis.line = element_line(colour = "black")) ## ----tidy=TRUE---------------------------------------------------------------- Selection<-VarSelection(Tara_Out,TaraOc,Crit="Rsquare",perc=0.95) ## ----tidy=TRUE---------------------------------------------------------------- library(reshape2) Ds<-as.data.frame(table(Variables(Selection))) ggplot(Ds, aes(x=Var1, y = Freq)) + geom_bar(stat="identity", fill="steelblue") + ggpubr::rotate_x_text(45) ## ----results='hide',eval=FALSE,tidy=TRUE-------------------------------------- # library(LinkHD) # data("Ruminotypes") ## ----results='hide', message=FALSE, warning=FALSE,tidy=TRUE---- library("LinkHD") library("ggplot2") options(width =50) ## ----tidy=TRUE---------------------------------- data("Ruminotypes") ## ----tidy=TRUE---------------------------------- Normalization<-lapply(list(Ruminotypes$`16_S`,Ruminotypes$Archaea,Ruminotypes$`18_S`) ,function(x){DataProcessing(x,Method="Compositional")}) Dataset<-Normalization names(Dataset)<-c("16_S","Archaea","18_S") ## ----tidy=TRUE---------------------------------- library(MultiAssayExperiment) library(LinkHD) data("Ruminotypes") Normalization<-lapply(list(Ruminotypes$`16_S`,Ruminotypes$Archaea,Ruminotypes$`18_S`), function(x){DataProcessing(x,Method="Compositional")}) Dataset<-Normalization names(Dataset)<-c("16_S","Archaea","18_S") Dataset<-MultiAssayExperiment(experiments = Dataset) ## ----tidy=TRUE---------------------------------- Output<-LinkData(Dataset,Distance=rep("euclidean",3),Scale = FALSE,Center=FALSE,nCluster = 3) #variance explained by compromise coordinates CompromisePlot(Output)+theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), panel.background = element_blank(), axis.line = element_line(colour="black")) ## ----eval=FALSE,tidy=TRUE----------------------- # B<-ExpressionSet() # Dataset1<-Normalization # A<-lapply(Dataset1,as.matrix) # B<-lapply(A,ExpressionSet) # Output_<-LinkData(B,Distance=rep("euclidean",3),Scale = FALSE,Center=FALSE,nCluster = 3) ## ----message=FALSE,tidy=TRUE-------------------- GlobalPlot(Output) ## ----tidy=TRUE---------------------------------- CorrelationPlot(Output) ## ----tidy=TRUE---------------------------------- mydata=compromise_coords(Output) #set colors cols <- c("1" = "red", "2" = "blue", "3" = "orange") mydata$fit.cluster<-as.factor(mydata$fit.cluster) p<-ggplot(mydata, aes(Dim.1,Dim.2)) + geom_point(aes(colour = fit.cluster)) p + scale_colour_manual(values = cols) + theme_classic() ## ----message=FALSE,tidy=TRUE-------------------- library( emmeans) Dat<-data.frame("y"=Ruminotypes$phenotype$ch4y,"E"=as.factor(mydata$fit.cluster)) s=lm(y~E,data=Dat) #summary(s) lsmeans(s, "E") pairs(lsmeans(s, "E")) ## ----message=FALSE,tidy=TRUE-------------------- Select_Var<-VarSelection(Output,Data=Dataset,Crit = "Rsquare",perc=0.9) ## ----message=FALSE,tidy=TRUE-------------------- M<-dAB(Output,Data=Dataset) ## ----message=FALSE,tidy=TRUE-------------------- #Select only the selected variables form Select_Var object SelectedB_16S<-Variables(Select_Var)[VarTable(Select_Var)=="16_S"] SelectedB_Archaea<-Variables(Select_Var)[VarTable(Select_Var)=="Archaea"] SelectedB_18S<-Variables(Select_Var)[VarTable(Select_Var)=="18_S"] ## ----message=FALSE,tidy=TRUE-------------------- Selec_16S<-as.matrix(Ruminotypes[[1]][,colnames(Dataset[[1]])%in%names(M[[1]])]) Selec_Archaea<-as.matrix(Ruminotypes[[2]][,colnames(Dataset[[2]])%in%names(M[[2]])]) Selec_18S<-as.matrix(Ruminotypes[[3]][,colnames(Dataset[[3]])%in%names(M[[3]])]) ## ----message=FALSE,tidy=TRUE-------------------- library(reshape2) SignTaxa<-OTU2Taxa(Selection=VarTable(Select_Var),TaxonInfo=Ruminotypes$Taxa_16S, tableName="16_S",AnalysisLev = "Family") melted_Table <- melt(SignTaxa$TotalUp1) ggplot(data =melted_Table, aes(x= reorder(rownames(melted_Table), -value), y=value)) + scale_fill_gradient(low = "steelblue", high = "#63D2E4",na.value="steelblue") + geom_bar(stat="identity", fill="#63D2E4")+ theme_classic()+ theme( axis.ticks = element_blank(),axis.title.x=element_blank(),axis.title.y = element_blank())+ ggpubr::rotate_x_text(90)+ theme(axis.text.x=element_text(size=rel(1.5))) ## ----message=FALSE,tidy=TRUE-------------------- SignTaxa<-OTU2Taxa(Selection=M,TaxonInfo=Ruminotypes$Taxa_16S,tableName="16_S",AnalysisLev = "Family") melted_Table <- melt(SignTaxa$TotalUp1) ggplot(data =melted_Table, aes(x= reorder(rownames(melted_Table), -value), y=value)) + scale_fill_gradient(low = "steelblue", high = "#FF0066",na.value="steelblue") + theme_classic()+ geom_bar(stat="identity", fill= "#56B4E9" )+ theme( axis.ticks = element_blank(),axis.title.x=element_blank(),axis.title.y = element_blank())+ ggpubr::rotate_x_text(90) + theme(axis.text.x=element_text(size=rel(0.9))) ## ----message=FALSE,tidy=TRUE-------------------- Rumin<-which(Ruminotypes$Taxa_16S[Ruminotypes$Taxa_16S$OTUID%in%SelectedB_16S,]$Family=="Ruminococcaceae") OTUs<-Ruminotypes$Taxa_16S[Ruminotypes$Taxa_16S$OTUID%in%SelectedB_16S,][Rumin,]$OTUID Rumin_Tot<-Dataset[['16_S']][,colnames(Dataset[['16_S']])%in%OTUs] RuminSum<-apply(Rumin_Tot,1,sum) Total<-data.frame(RuminSum,"Group"=as.factor(compromise_coords(Output)$fit.cluster)) p <- ggplot(Total, aes(x=Group, y=RuminSum, fill=Group)) + geom_boxplot()+ scale_fill_manual(values=c("Tomato", "Orange", "DodgerBlue"))+ theme_classic()+ ggtitle(paste("Differences at level",rownames(melted_Table)[16])) p ## ----message=FALSE,tidy=TRUE-------------------- Christ<-which(Ruminotypes$Taxa_16S[Ruminotypes$Taxa_16S$OTUID%in%SelectedB_16S,]$Family==rownames(melted_Table)[14]) OTUs<-Ruminotypes$Taxa_16S[Ruminotypes$Taxa_16S$OTUID%in%SelectedB_16S,][Christ,]$OTUID Christ_Tot<-Dataset[['16_S']][,colnames(Dataset[['16_S']])%in%OTUs] ChristSum<-apply(Christ_Tot,1,sum) Total<-data.frame(ChristSum,"Group"=as.factor(compromise_coords(Output)$fit.cluster)) p <- ggplot(Total, aes(x=Group, y=ChristSum, fill=Group)) + geom_boxplot()+ scale_fill_manual(values=c("Tomato", "Orange", "DodgerBlue"))+ theme_classic()+ ggtitle(paste("Differences at level",rownames(melted_Table)[14])) p ## ----message=FALSE,tidy=TRUE-------------------- ggplot(Total, aes(x=Group, y=rownames(Total), fill=ChristSum)) + geom_tile() + scale_fill_gradient(low = "steelblue", high = "#FF0066",na.value="steelblue") + ggpubr::rotate_x_text(45)+ theme_classic() + theme( axis.ticks = element_blank(),axis.title.x=element_blank(),axis.title.y = element_blank())+ theme(axis.text.y=element_text(size=rel(0.6))) ## ----message=FALSE,tidy=TRUE-------------------- Succi<-which(Ruminotypes$Taxa_16S[Ruminotypes$Taxa_16S$OTUID%in%SelectedB_16S,]$Family=="Succinivibrionaceae") OTUs<-Ruminotypes$Taxa_16S[Ruminotypes$Taxa_16S$OTUID%in%SelectedB_16S,][Succi,]$OTUID Succi_Tot<-Dataset[['16_S']][,colnames(Dataset[['16_S']])%in%OTUs] SucciSum<-apply(Succi_Tot,1,sum) Total<-data.frame(SucciSum,"Group"=as.factor(compromise_coords(Output)$fit.cluster)) Total<-Total[- which(Total[,1]>1),] p <- ggplot(Total, aes(x=Group, y=SucciSum, fill=Group)) + geom_boxplot()+ scale_fill_manual(values=c("#e46563", "#357fe3", "#75e463"))+ theme_classic()+ ggtitle(paste("Differences at level","Succinivibrionaceae"))+ theme(axis.text.x=element_text(size=rel(1.5)),legend.title=element_text(size=rel(1.5)), legend.text=element_text(size=rel(1.5))) p ## ----message=FALSE,tidy=TRUE-------------------- SignTaxa<-OTU2Taxa(Selection=VarTable(Select_Var), TaxonInfo=Ruminotypes$Taxa_18S,tableName="18_S") melted_Table <- melt(SignTaxa$TotalUp1) ggplot(data =melted_Table, aes(x= reorder(rownames(melted_Table), -value), y=value)) + scale_fill_gradient(low = "steelblue", high = "#B848DB",na.value="steelblue") + geom_bar(stat="identity", fill="#B848DB")+ theme_classic()+ theme( axis.ticks = element_blank(),axis.title.x=element_blank(),axis.title.y = element_blank())+ ggpubr::rotate_x_text(90)+ theme(axis.text.x=element_text(size=rel(1.5))) ## ----message=FALSE,tidy=TRUE-------------------- SignTaxa<-OTU2Taxa(Selection=M,TaxonInfo=Ruminotypes$Taxa_18S,tableName="18_S") melted_Table <- melt(SignTaxa$TotalUp1[SignTaxa$TotalUp1>1]) ggplot(data =melted_Table, aes(x= reorder(rownames(melted_Table), -value), y=value)) + scale_fill_gradient(low = "steelblue", high = "#FF0066",na.value="steelblue") + geom_bar(stat="identity", fill= "#56B4E9" )+ theme_classic()+ theme( axis.ticks = element_blank(),axis.title.x=element_blank(),axis.title.y = element_blank())+ ggpubr::rotate_x_text(90)+ theme(axis.text.x=element_text(size=rel(0.8))) ## ----message=FALSE,tidy=TRUE-------------------- Dip<-which(Ruminotypes$Taxa_18S[Ruminotypes$Taxa_18S$OTUID%in%SelectedB_18S,]$Genus=="Diplodinium") OTUs<-Ruminotypes$Taxa_18S[Ruminotypes$Taxa_18S$OTUID%in%SelectedB_18S,][Dip,]$OTUID Dip_Tot<-Dataset[['18_S']][,colnames(Dataset[['18_S']])%in%OTUs] DipSum<-apply(Dip_Tot,1,sum) Total<-data.frame(DipSum,"Group"=as.factor(compromise_coords(Output)$fit.cluster)) p <- ggplot(Total, aes(x=Group, y=DipSum, fill=Group)) + geom_boxplot()+ scale_fill_manual(values=c("Tomato", "Orange", "DodgerBlue"))+ theme_classic()+ ggtitle(paste("Differences at level","Diplodinium")) p ## ----message=FALSE,tidy=TRUE-------------------- ggplot(Total, aes(x=Group, y=rownames(Total), fill=DipSum)) + geom_tile() + scale_fill_gradient(low = "steelblue", high = "#FF0066",na.value="steelblue") + ggpubr::rotate_x_text(45)+ theme_classic() + theme( axis.ticks = element_blank(),axis.title.x=element_blank(),axis.title.y = element_blank())+ theme(axis.text.y=element_text(size=rel(0.6))) ## ----message=FALSE,warning=FALSE,tidy=TRUE------ SignTaxa<-OTU2Taxa(Selection=VarTable(Select_Var),TaxonInfo=Ruminotypes$Taxa_Archaea,tableName="Archaea") melted_Table <- melt(SignTaxa$TotalUp1) ggplot(data =melted_Table, aes(x= reorder(rownames(melted_Table), -value), y=value)) + scale_fill_gradient(low = "steelblue", high = "#46D6AD",na.value="steelblue") + geom_bar(stat="identity", fill="#46D6AD")+ theme_classic()+ theme( axis.ticks = element_blank(),axis.title.x=element_blank(),axis.title.y = element_blank())+ ggpubr::rotate_x_text(90)+ theme(axis.text.x=element_text(size=rel(1.5))) ## ----message=FALSE,warning=FALSE,tidy=TRUE------ SignTaxa<-OTU2Taxa(Selection=M,TaxonInfo=Ruminotypes$Taxa_Archaea,tableName="Archaea") melted_Table <- melt(SignTaxa$TotalUp1) ggplot(data =melted_Table, aes(x= reorder(rownames(melted_Table), -value), y=value)) + scale_fill_gradient(low = "steelblue", high = "#FF0066",na.value="steelblue") + geom_bar(stat="identity", fill= "#56B4E9" )+ theme_classic()+ theme( axis.ticks = element_blank(),axis.title.x=element_blank(),axis.title.y = element_blank())+ ggpubr::rotate_x_text(90)+ theme(axis.text.x=element_text(size=rel(0.9))) ## ----message=FALSE,warning=FALSE,tidy=TRUE------ Met<-which(Ruminotypes$Taxa_Archaea[Ruminotypes$Taxa_Archaea$OTU%in%SelectedB_Archaea,]$Genus==levels(Ruminotypes$Taxa_Archaea$Genus)[4]) OTUs<-Ruminotypes$Taxa_Archaea[Ruminotypes$Taxa_Archaea$OTU%in%SelectedB_Archaea,][Met,]$OTU Met_Tot<-Dataset[['Archaea']][,colnames(Dataset[['Archaea']])%in%OTUs] MetSum<-apply(Met_Tot,1,sum) Total<-data.frame(MetSum,"Group"=as.factor(compromise_coords(Output)$fit.cluster)) p <- ggplot(Total, aes(x=Group, y=MetSum, fill=Group)) + geom_boxplot()+ scale_fill_manual(values=c("Tomato", "Orange", "DodgerBlue"))+ theme_classic()+ ggtitle(paste("Differences at level","Methanobrevibacter")) p ## ----message=FALSE,warning=FALSE,tidy=TRUE------ ggplot(Total, aes(x=Group, y=rownames(Total), fill=MetSum)) + geom_tile() + scale_fill_gradient(low = "steelblue", high = "#FF0066",na.value="steelblue") + ggpubr::rotate_x_text(45)+ theme_classic() + theme( axis.ticks = element_blank(),axis.title.x=element_blank(),axis.title.y = element_blank())+ theme(axis.text.y=element_text(size=rel(0.6))) ## ----message=FALSE,warning=FALSE,eval=FALSE,tidy=TRUE---- # #library("readxl") # #library("XLConnect") # #library("openxlsx") # # # Download the list from supplementary material: # #temp = tempfile(fileext = ".xls") # #dataURL <- "https://onlinelibrary.wiley.com/action/downloadSupplement?doi=10.1111%2Fjbg.12427&file=jbg12427-sup-0004-TableS4.xls" # #download.file(dataURL, destfile=temp, mode='wb') # #SPLS_DA<-read_excel(temp,sheet=1) # #SPLS<-read_excel(temp,sheet=2) # # #Create a list with selected OTU. Note that if you have only one table, the list should have lenght equal to 1. # # #Selected_OTU<-list("SPLS"=SPLS$OTU,"SPLS_DA"=SPLS_DA$OTU) # # OTUs have to be the lists names: # #names(Selected_OTU$SPLS)<-SPLS$OTU # #names(Selected_OTU$SPLS_DA)<-SPLS_DA$OTU # # # #TaxonInfo is the same than 16_S taxon info into Ruminotypes dataset: # #Apply OTU2Taxa # #SignTaxa<-OTU2Taxa(Selection=Selected_OTU,TaxonInfo=Ruminotypes$Taxa_16S,tableName="SPLS",AnalysisLev = "Family") # #SignTaxa2<-OTU2Taxa(Selection=Selected_OTU,TaxonInfo=Ruminotypes$Taxa_16S,tableName="SPLS_DA",AnalysisLev = "Family") # # # # #melted_Table <- melt(SignTaxa2$TotalUp1) # #ggplot(data =melted_Table, aes(x= reorder(rownames(melted_Table), -value), y=value)) + # #scale_fill_gradient(low = "steelblue",high = "#FF0066",na.value="steelblue") + # # theme_classic()+ # # geom_bar(stat="identity", fill= "#ff661a" )+ # # theme( axis.ticks = element_blank(),axis.title.x=element_blank(),axis.title.y = element_blank())+ # # ggpubr::rotate_x_text(90) + # # theme(axis.text.x=element_text(size=rel(0.9)))