3.1 Loading data from ncbi GEO

Load the files and make a list of the matrix:

When reading the normalized expression values the function check whether there are rows that their median absolute deviation (MAD) value equal to zero and remove these rows. This step is important in order to continue to the next step of indexing the data.

library(GEOquery)

gse2077 <- getGEO('GSE2077')
gseAssays <- Biobase::assayData(gse2077)
gseExprs <- Biobase::assayDataElement(gseAssays[[1]][,c(1:6)],'exprs')

# list of the time series tables use only 6 samples
gseList <- lapply(1:6,function(x){data.frame(Genes = rownames(gseExprs),Value = gseExprs[,x])})
names(gseList) <- colnames(gseExprs)

Build the expression matrix from the list of matrices:

rts <- readTSGE(gseList,labels = c("0h","6h","12h","24h","48h","72h")) 

3.2 Loading data from files

Here is an example of how to load the time-series files from your directory; data were imported from the ctsGE package:

data_dir <- system.file("extdata", package = "ctsGE")
files <- dir(path=data_dir,pattern = "\\.xls$")

Building from a directory:

rts <- readTSGE(files, path = data_dir, labels = c("0h","6h","12h","24h","48h","72h") )

ctsGE Object summary:

names(rts)

## [1] "tsTable"    "samples"    "tags"       "timePoints"

rts$timePoints

## [1] 6

head(rts$samples)

## [1] "0h"  "6h"  "12h" "24h" "48h" "72h"

head(rts$tags)

## [1] "1000_at"   "1001_at"   "1002_f_at" "1003_s_at" "1004_at"   "1005_at"

head(rts$tsTable)

4.1 Checking for optimal cutoff

The function PreparingTheIndexes generates an expression index (i.e., a sequence of 1,-1, and 0) that represents the expression pattern along time points for each gene. Setting different limits for the center era (with the parameter cutoff) will change the index for each gene-expression profile and consequently, the number of genes in each index group. Following the idea that instead of filtering “irrelevant” genes to reduce the noise, the clustering will be performed on small gene groups, one would like to choose a cutoff value, that will minimize the number of genes in each group, i.e., generate index groups of equal size. The test for equality is performed by calculating the chi-squared value from a comparison between the number of genes in the index groups and the null hypothesis that all index groups are equal. The test is performed for all the cutoffs in the range and the cutoff that gives the minimal chi-squared value is the most likely to generate equal index groups. The range of the cutoff values is given by min_cutoff and max_cutoff arguments. However, by setting the same value to min and max parameters one can define a cutoff regardless of what was suggested by the function.

prts <- PreparingTheIndexes(x = rts, min_cutoff=0.5, max_cutoff=0.7, mad.scale = TRUE)

cutoff = 0.55 is the optimal cutoff with the lowest chi-squared value

prts$cutoff

## [1] 0.55

With cutoff =
0.55, most of the genes were assigned to indexes with three or two zeros,
indicating a variety of expression patterns.

4.2 Preparing the indexes for the data

A cutoff = 0.55 was chosen

prts <- PreparingTheIndexes(x = rts, mad.scale = TRUE)
names(prts)

## [1] "tsTable"    "samples"    "tags"       "timePoints" "scaled"
## [6] "index"      "cutoff"

Gene expression after standardization: head(prts$scaled)

Gene expression indexing with cutoff = 0.55: head(prts$index)

4.3 Clustering each index with K-means

The clustering is done with K-means. To choose an optimal k for K-means clustering, the Elbow method was applied⁴, this method looks at the percentage of variance explained as a function of the number of clusters: the chosen number of clusters should be such that adding another cluster does not give much better modeling of the data. First, the ratio of the within-cluster sum of squares (WSS) to the total sum of squares (TSS) is computed for different values of k (i.e., 1, 2, 3 …). The WSS, also known as sum of squared error (SSE), decreases as k gets larger. The Elbow method chooses the k at which the SSE decreases abruptly. This happens when the computed value of the WSS-to-TSS ratio first drops from 0.2.

\(\frac{WSS}{TSS} < 0.2\)

ClustIndexes <- ClustIndexes(prts, scaling = TRUE)
names(ClustIndexes)
# table of the index and the recommended k that were found by the function
head(ClustIndexes$optimalK)

# Table of clusters index for each gene
head(ClustIndexes$ClusteredIdxTable) 

Running kmeans and calculating the optimal k for each one of the indexes in the data could take a long time. To shorten the procedure the user can skip this step altogether and directly view a specific index and its clusters by running either the PlotIndexesClust() or the ctsGEShinyApp() function

5.1 Export gene expression table to file

Genes expression data in index ‘1100-1-1’ separated into clsuter can be exported with R function write.table:

write.table(indexPlot[[1]], file, sep = "\t")

6.1 Screenshots of the GUI

library(shiny)
library(DT)
ctsGEShinyApp(rts)