In this vignette, we describe how to use the NewmanOmics Paired and
Banked tests to analyze gene expression data from a single sample.
Getting Started
As usual, we start by loading the package:
The package contains paired tumor and normal samples from patients
with head and neck cancer. these came from a study that was submitted to
the Gene Expression Omnibus.
data(GSE6631)
dim(GSE6631)
## [1] 2000 44
## Normal.mucosa.1 Cancer.1 Normal.mucosa.2 Cancer.2
## 34155_s_at 26.42586 22.19725 22.13673 18.66223
## 34281_at 334.29232 382.92879 393.40014 509.30754
## 39125_at 258.62695 290.06060 268.97994 220.16837
## 37276_at 45.65556 38.86692 34.77368 33.40627
## 1519_at 423.26690 366.40731 308.62338 550.34888
As we can see, this consists of (normalized) Affymetrix microarray
data. The odd numbered columns are derived from normal mucosa, and the
even numbered columns are derived from paired tumor samples.
Before proceeding, we are going to log-transform the data.
HN <- log2(1 + GSE6631)
boxplot(HN, col=c("forestgreen", "dodgerblue"))
The figure suggests that the
the data have been reasonably normalized, and that it is unlikely to be
overwhelmed by artifacts.
Paired Statistic
To illustrate the Newman Paired test, we are going to use only one
sample.
HN1 <- HN[, 1:2]
result1 <- pairedStat(HN1, pairing = c(-1,1))
summary(result1@nu.statistics)
## Cancer.1
## Min. :5.842e-04
## 1st Qu.:4.158e-01
## Median :9.885e-01
## Mean :1.417e+00
## 3rd Qu.:1.835e+00
## Max. :1.744e+01
summary(result1@p.values)
## Cancer.1
## Min. :0.0000
## 1st Qu.:0.3005
## Median :0.5777
## Mean :0.5464
## 3rd Qu.:0.8150
## Max. :0.9995
We can create a histogram of the per-gene (empirical) p-values
We can also produce an
“M-versus-A” plot of the data.
Alternate Inputs
The pairedStat function has flexible inputs, allowing you to
store the data in various ways. Here we run the algorithm for three
pairs, with an explicit pairing vector.
result2 <- pairedStat(HN[, 1:6], pairing=c(-1, 1, -2, 2, -3, 3))
summary(result2@nu.statistics)
## Cancer.1 Cancer.2 Cancer.3
## Min. :5.842e-04 Min. :9.836e-04 Min. : 0.001682
## 1st Qu.:4.158e-01 1st Qu.:4.591e-01 1st Qu.: 0.369775
## Median :9.885e-01 Median :1.021e+00 Median : 0.878206
## Mean :1.417e+00 Mean :1.417e+00 Mean : 1.419570
## 3rd Qu.:1.835e+00 3rd Qu.:1.780e+00 3rd Qu.: 1.765563
## Max. :1.744e+01 Max. :1.526e+01 Max. :16.572431
summary(result2@p.values)
## Cancer.1 Cancer.2 Cancer.3
## Min. :0.0000 Min. :0.0000 Min. :0.0000
## 1st Qu.:0.3002 1st Qu.:0.3153 1st Qu.:0.3191
## Median :0.5764 Median :0.5651 Median :0.6185
## Mean :0.5463 Mean :0.5376 Mean :0.5629
## 3rd Qu.:0.8151 3rd Qu.:0.7961 3rd Qu.:0.8350
## Max. :0.9988 Max. :0.9997 Max. :0.9995
We can also input the same data as a pair of matrices.
normals <- HN[, c(1,3,5)]
tumors <- HN[, c(2,4,6)]
result3 <- pairedStat(normals, tumors)
summary(result3@nu.statistics)
## Cancer.1 Cancer.2 Cancer.3
## Min. :5.842e-04 Min. :9.836e-04 Min. : 0.001682
## 1st Qu.:4.158e-01 1st Qu.:4.591e-01 1st Qu.: 0.369775
## Median :9.885e-01 Median :1.021e+00 Median : 0.878206
## Mean :1.417e+00 Mean :1.417e+00 Mean : 1.419570
## 3rd Qu.:1.835e+00 3rd Qu.:1.780e+00 3rd Qu.: 1.765563
## Max. :1.744e+01 Max. :1.526e+01 Max. :16.572431
summary(result3@p.values)
## Cancer.1 Cancer.2 Cancer.3
## Min. :0.0000 Min. :0.0000 Min. :0.0000
## 1st Qu.:0.3009 1st Qu.:0.3158 1st Qu.:0.3199
## Median :0.5762 Median :0.5640 Median :0.6189
## Mean :0.5461 Mean :0.5373 Mean :0.5628
## 3rd Qu.:0.8144 3rd Qu.:0.7954 3rd Qu.:0.8349
## Max. :0.9999 Max. :1.0000 Max. :0.9997
Or we can input the same data as a list of paired samples.
listOfPairs <- list(HN[,1:2], HN[,3:4], HN[,5:6])
result4 <- pairedStat(listOfPairs)
summary(result4@nu.statistics)
## Cancer.1 Cancer.2 Cancer.3
## Min. :5.842e-04 Min. :9.836e-04 Min. : 0.001682
## 1st Qu.:4.158e-01 1st Qu.:4.591e-01 1st Qu.: 0.369775
## Median :9.885e-01 Median :1.021e+00 Median : 0.878206
## Mean :1.417e+00 Mean :1.417e+00 Mean : 1.419570
## 3rd Qu.:1.835e+00 3rd Qu.:1.780e+00 3rd Qu.: 1.765563
## Max. :1.744e+01 Max. :1.526e+01 Max. :16.572431
summary(result4@p.values)
## Cancer.1 Cancer.2 Cancer.3
## Min. :0.0000 Min. :0.0000 Min. :0.0000
## 1st Qu.:0.3013 1st Qu.:0.3154 1st Qu.:0.3199
## Median :0.5780 Median :0.5639 Median :0.6206
## Mean :0.5463 Mean :0.5376 Mean :0.5632
## 3rd Qu.:0.8164 3rd Qu.:0.7949 3rd Qu.:0.8360
## Max. :0.9996 Max. :0.9990 Max. :0.9996
Banked Statistic
A completely different approach to personalized transcriptomics is to
compare individual samples to a “bank” of known normals.
normals <- HN[, seq(1, ncol(HN), 2)] # odds are normal
tumors <- HN[, seq(2, ncol(HN), 2)] # evens are tumor
bank <- createBank(normals)
result5 <- bankStat(bank, tumors[,1,drop=FALSE])
summary(result5$nu.statistics)
## Cancer.1
## Min. :-9.6255
## 1st Qu.:-0.5540
## Median : 0.2797
## Mean : 0.1027
## 3rd Qu.: 0.9911
## Max. : 6.7524
summary(result5$p.values)
## Cancer.1
## Min. :0.0000
## 1st Qu.:0.2898
## Median :0.6101
## Mean :0.5569
## 3rd Qu.:0.8392
## Max. :1.0000
hist(result5$p.values, breaks=101)
