Nothing
## ---- echo = FALSE, message = FALSE-----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
library(myTAI)
options(width = 750)
knitr::opts_chunk$set(
comment = "#>",
error = FALSE,
tidy = FALSE)
## ----eval=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# library(myTAI)
#
# # load an example PhyloExpressionSet stored in the myTAI package
# data(PhyloExpressionSetExample)
#
# # look at the standardized data set format
# head(PhyloExpressionSetExample, 3)
## ----eval=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # TAI profile of correctly assigned Phylostrata
# TAI(PhyloExpressionSetExample)
## ---- fig.width= 7, fig.height= 5,eval=FALSE--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# data(PhyloExpressionSetExample)
# # Visualize the TAI profile of correctly assigned Phylostrata
# PlotSignature( ExpressionSet = PhyloExpressionSetExample,
# p.value = FALSE )
## ----eval=FALSE, echo=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # randomly permute the Phylostratum assignment of all genes
# randomPhyloExpressionSetExample <- data.frame( sample(PhyloExpressionSetExample[ , "Phylostratum"]),
# PhyloExpressionSetExample[, 2:9] )
#
# # TAI profile based on randomly assigned Phylostrata
# TAI(randomPhyloExpressionSetExample)
## ---- fig.width= 7, fig.height= 5,eval=FALSE--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # Visualize the TAI profile based on randomly assigned Phylostrata
# PlotSignature( ExpressionSet = randomPhyloExpressionSetExample,
# p.value = FALSE )
## ----eval=FALSE, echo=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # Variance of the TAI profile based on correctly assigned Phylostrata
# var(TAI(PhyloExpressionSetExample))
## ----eval=FALSE, echo=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # Variance of the TAI profile based on randomly assigned Phylostrata
# var(TAI(randomPhyloExpressionSetExample))
## ----eval=FALSE, echo=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # Compute 1000 TAI profiles based on randomly assigned Phylostrata
# randomTAIs <- bootMatrix( ExpressionSet = PhyloExpressionSetExample,
# permutations = 1000 )
#
# head(randomTAIs)
## ---- fig.width= 7, fig.height= 5,eval=FALSE--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # compute the variance of the random TAI profile for each row
# variance_vector <- apply(randomTAIs, 1 , var)
#
# # and visualize the distribution of variances
# hist(variance_vector, breaks = 100)
## ---- fig.width= 7, fig.height= 5,eval=FALSE--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # variance of the TAI profile based on correctly assigned Phylostrata
# var_real <- var(TAI(PhyloExpressionSetExample))
#
# # visualize the distribution of variances
# hist( x = c(variance_vector,var_real),
# breaks = 100,
# xlab = "variance",
# main = "Histogram of variance_vector" )
#
# # and plot a red line at the position where we can find the
# # real variance
# abline(v = var_real, lwd = 5, col = "red")
## ---- fig.width= 7, fig.height= 5,eval=FALSE--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # install.packages("fitdistrplus")
#
# # plot a Cullen and Frey graph
# fitdistrplus::descdist(variance_vector)
## ---- fig.width= 7, fig.height= 5,eval=FALSE--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # estimate the parameters: shape and rate using 'moment matching estimation'
# gamma_MME <- fitdistrplus::fitdist(variance_vector,distr = "gamma", method = "mme")
# # estimate shape:
# shape <- gamma_MME$estimate[1]
# # estimate the rate:
# rate <- gamma_MME$estimate[2]
#
# # define an expression written as function as input for the curve() function
# gamma_distr <- function(x){ return(dgamma(x = x, shape = shape, rate = rate)) }
#
# # plot the density function and the histogram of variance_vector
# curve( expr = gamma_distr,
# xlim = c(min(variance_vector),max(c(variance_vector,var_real))),
# col = "steelblue",
# lwd = 5,
# xlab = "Variances",
# ylab = "Frequency" )
#
# # plot the histogram of variance_vector
# histogram <- hist(variance_vector,prob = TRUE,add = TRUE, breaks = 100)
# rug(variance_vector)
#
# # plot a red line at the position where we can find the real variance
# abline(v = var_real, lwd = 5, col = "red")
## ----eval=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # p-value of var_real
# pgamma(var_real, shape = shape,rate = rate, lower.tail = FALSE)
## ----eval=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # Perform the FlatLineTest
# FlatLineTest( ExpressionSet = PhyloExpressionSetExample,
# permutations = 1000 )
## ---- fig.width= 7, fig.height= 3,eval=FALSE--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
#
# # perform the FlatLineTest and investigate the goodness of the test statistic
# FlatLineTest( ExpressionSet = PhyloExpressionSetExample,
# permutations = 1000,
# plotHistogram = TRUE )
#
## ----eval=FALSE, echo=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# library(myTAI)
#
# # load an example PhyloExpressionSet stored in the myTAI package
# data(PhyloExpressionSetExample)
#
# # look at the standardized data set format
# head(PhyloExpressionSetExample, 3)
## ----eval=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # TAI profile of correctly assigned Phylostrata
# TAI(PhyloExpressionSetExample)
## ---- fig.width= 7, fig.height= 5,eval=FALSE--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # visualize the TAI profile of correctly assigned Phylostrata
# PlotSignature( ExpressionSet = PhyloExpressionSetExample,
# p.value = FALSE )
## ---- fig.width= 7, fig.height= 5,eval=FALSE--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # Perform the Reductive Hourglass Test
# ReductiveHourglassTest( ExpressionSet = PhyloExpressionSetExample,
# modules = list(early = 1:2, mid = 3:5, late = 6:7),
# lillie.test = TRUE )
#
## ---- fig.width= 7, fig.height= 6,eval=FALSE--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # perform the Reductive Hourglass Test and plot the test statistic
# ReductiveHourglassTest( ExpressionSet = PhyloExpressionSetExample,
# modules = list(early = 1:2, mid = 3:5, late = 6:7),
# plotHistogram = TRUE,
# lillie.test = TRUE )
#
## ----eval=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# library(myTAI)
#
# # load an example PhyloExpressionSet stored in the myTAI package
# data(PhyloExpressionSetExample)
#
# # look at the standardized data set format
# head(PhyloExpressionSetExample, 3)
## ----eval=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # TAI profile of correctly assigned Phylostrata
# TAI(PhyloExpressionSetExample)
## ---- fig.width= 7, fig.height= 5,eval=FALSE--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # Visualize the TAI profile of correctly assigned Phylostrata
# PlotSignature( ExpressionSet = PhyloExpressionSetExample,
# p.value = FALSE )
#
## ---- fig.width= 7, fig.height= 5,eval=FALSE--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
#
# # Perform the Reductive Early Conservation Test
# EarlyConservationTest( ExpressionSet = PhyloExpressionSetExample,
# modules = list(early = 1:2, mid = 3:5, late = 6:7),
# lillie.test = TRUE )
#
## ---- fig.width= 7, fig.height= 6,eval=FALSE--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # perform the Reductive Early Conservation Test and plot the test statistic
# EarlyConservationTest( ExpressionSet = PhyloExpressionSetExample,
# modules = list(early = 1:2, mid = 3:5, late = 6:7),
# plotHistogram = TRUE,
# lillie.test = TRUE )
#
## ----eval=FALSE, echo=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# library(myTAI)
#
# # load an example PhyloExpressionSet stored in the myTAI package
# data(PhyloExpressionSetExample)
#
# # look at the standardized data set format
# head(PhyloExpressionSetExample, 3)
## ----eval=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # TAI profile of correctly assigned Phylostrata
# TAI(PhyloExpressionSetExample)
## ---- fig.width= 7, fig.height= 5,eval=FALSE--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # visualize the TAI profile of correctly assigned Phylostrata
# PlotSignature( ExpressionSet = PhyloExpressionSetExample,
# p.value = FALSE )
## ---- fig.width= 7, fig.height= 5,eval=FALSE--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # Perform the Reverse Hourglass Test
# ReverseHourglassTest( ExpressionSet = PhyloExpressionSetExample,
# modules = list(early = 1:2, mid = 3:5, late = 6:7),
# lillie.test = TRUE )
#
## ---- fig.width= 7, fig.height= 6,eval=FALSE--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # perform the Reverse Hourglass Test and plot the test statistic
# ReverseHourglassTest( ExpressionSet = PhyloExpressionSetExample,
# modules = list(early = 1:2, mid = 3:5, late = 6:7),
# plotHistogram = TRUE,
# lillie.test = TRUE )
#
## ---- fig.width= 7, fig.height= 5-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
library(myTAI)
data(PhyloExpressionSetExample)
# a simple example is to transform the gene expression levels of a given PhyloExpressionSet
# using a sqrt or log2 transformation
PES.sqrt <- tf(PhyloExpressionSetExample, sqrt)
head(PES.sqrt)
## ---- fig.width= 7, fig.height= 5, eval=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
#
# PES.log2 <- tf(PhyloExpressionSetExample, log2)
#
# head(PES.log2)
## ---- fig.width= 7, fig.height= 5,eval=FALSE--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# # in case a given PhyloExpressionSet already stores gene expression levels
# # that are log2 transformed and need to be re-transformed to absolute
# # expression levels, to perform subsequent phylotranscriptomics analyses
# # (that are defined for absolute expression levels),
# # one can re-transform a PhyloExpressionSet like this:
#
# PES.absolute <- tf(PES.log2 , function(x) 2^x)
#
# # which should be the same as PhyloExpressionSetExample :
# head(PhyloExpressionSetExample)
# head(PES.absolute)
## ---- fig.width= 7, fig.height= 5, eval=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# data(PhyloExpressionSetExample)
# # plotting the TAI using log2 transformed expression levels
# # and performing the Flat Line Test to obtain the p-value
# PlotSignature( ExpressionSet = tf(PhyloExpressionSetExample, log2),
# TestStatistic = "FlatLineTest",
# xlab = "Ontogeny",
# ylab = "TAI" )
#
#
## ---- fig.width= 7, fig.height= 5, eval=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# data(PhyloExpressionSetExample)
# # plotting the TAI using sqrt transformed expression levels
# # and performing the Flat Line Test to obtain the p-value
# PlotSignature( ExpressionSet = tf(PhyloExpressionSetExample, sqrt),
# TestStatistic = "FlatLineTest",
# xlab = "Ontogeny",
# ylab = "TAI" )
## ---- fig.width= 7, fig.height= 5, eval=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# data(DivergenceExpressionSetExample)
# # plotting the TDI using log2 transformed expression levels
# # and performing the Flat Line Test to obtain the p-value
# PlotSignature( ExpressionSet = tf(DivergenceExpressionSetExample, log2),
# TestStatistic = "FlatLineTest",
# xlab = "Ontogeny",
# ylab = "TDI" )
## ---- fig.width= 7, fig.height= 5, eval=FALSE-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
# data(DivergenceExpressionSetExample)
# # plotting the TDI using sqrt transformed expression levels
# # and performing the Flat Line Test to obtain the p-value
# PlotSignature( ExpressionSet = tf(DivergenceExpressionSetExample, sqrt),
# TestStatistic = "FlatLineTest",
# xlab = "Ontogeny",
# ylab = "TDI" )
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.