Transcript exptbio - Pages - University of Wisconsin

A Graphical Investigation
of Some Microarray Experiments
Brian S. Yandell
Statistics , Horticulture, Biometry,
University of Wisconsin-Madison
www.stat.wisc.edu/~yandell/statgen
© BS Yandell
Experimental Biology 2002
1
Key Questions
• Why design microarray experiments? (Kerr Churchill)
– chips and samples are expensive
– design experiment for one or a few genes (want true replication)
• Are typical statistical assumptions warranted?
– how to transform to symmetry (near normal)?
– how does the variance change? by gene? with abundance?
• How do we combine data analysis across multiple genes?
– differential expression pattern changes with abundance
• how to keep potentially important low abundance genes?
– noise pattern changes with abundance
• How can we map gene expression?
– use pattern of expression as one or more quantitative traits
• Illustrate ideas with experiments from Attie Lab
© BS Yandell
Experimental Biology 2002
2
But what about MY technology?
• talk focuses on Affymetrix mouse chips
– 13,000+ mRNA fragments (11,000+ genes)
–  = mean(PM) - mean(MM) adjusted expression levels
• adaptable to other molecular data types
– genome, proteome, metabolome, megagenome, virome
• adaptable to emerging “micro-array” technologies
–
–
–
–
spotted arrays (Brown Botstein 1999)
micro-beads (www.lynxgen.com)
surface plasmon resonance (Nelson et al. Corn 2001)
maskless array synthesizer (www.nimblegen.com)
© BS Yandell
Experimental Biology 2002
3
Design: learning by experience
B6
• fat and obesity
F1 BTBR
lean
– lean vs. obese
– 3 “strains”
– no replicates: 4 mice per chip
• SCD knockout mouse
obese
knockout
– 5 replicates: 1 mouse per chip
– knockout vs. wild type
– 8 error degrees of freedom
wild type
• fat, liver, muscle, islet tissues
– 2 strains, 4 tissues
– 2 replicates: 2 mice per chip
– 8 error degrees of freedom
B6
B6
fat liver muscle islet
BTBR
BTBR
© BS Yandell
Experimental Biology 2002
4
How are obesity & diabetes related?
• focus on adipose (fat) tissue
– whole-body fuel partitioning
– Nadler et al. (2000) PNAS
B6
• 6 conditions in 2x3 factorial
– lean vs. obese
– strains B6, BTBR, F1 cross
• pseudo-replication = subsampling
– only 1 chip per condition
– 4 mice pooled per chip
F1
BTBR
lean
obese
• increase precision per chip
• but reduce power to detect change
• combine data across genes
– no way to infer differences otherwise
– noise decreases with average intensity
© BS Yandell
Experimental Biology 2002
5
SCD knockout experiment
• single gene knockout
knockout
– stearoyl-CoA desaturase-1
• experimental design
wild type
– knockout vs. wild type mice
– 5 mice per group, 1 chip per mouse
– dChip recalc of  = PM-MM
• have gene-specific replication
– estimate noise from replicates within groups
• compare genes in functional groups
– up or down regulation?
© BS Yandell
Experimental Biology 2002
6
Diabetes action in whole body
• tissues important for diabetes
– fat, muscle, liver, islets
– focus on fat & liver here
B6
B6
• two obese strains
fat
liver
muscle
islet
– BTBR diabetic
– B6 non-diabetic
BTBR
• experimental design
–
–
–
–
BTBR
only 16 Affymetrix chips
2 replicates each tissue*genotype condition
4 mice per condition in pools of 2 per chip
some benefits of pooling & independent replication
© BS Yandell
Experimental Biology 2002
7
How do we infer strain differences?
1.3
1.1
1.6
1.4
WT
KO
WT
SCD1
1
2
1.85
PPARalpha
-1
0
1.75
1.65
ACO
KO
PEPCK
-1.0
KO
WT
KO
WT
G6Pase
2.10
0.0
1.50 1.60 1.70
WT
2.40
WT
2.25
KO
© BS Yandell
KO
PPARgamma
knockout
wild type
GPAT
1.8
-0.5
-0.8
-1.1
WT
1.0
high noise
signal unreliable
test inconclusive
KO
-1.0
noise negligible?
-1.2
same mean
flat line
FAS
-1.4
strain differences?
SREBP1
KO
Experimental Biology 2002
WT
KO
WT
8
Why is noise so important?
• is differential expression “signal” large relative to “noise”?
– signal is difference across conditions of interest
• lean vs. obese, knockout vs. wild, B6 vs. BTBR
– noise assessed by “true” replicates, not pseudo-replicates
• sources of noise
– conditions: mouse, strain, tissue
• can vary with mRNA abundance, gene-specific features
– materials: chip, mRNA sample preparation
• hybridization and reading mechanics
– watch out for pseudo-replication
• pooled mRNA from multiple mice on one chip
• multiple chips from same mRNA source
• experimental unit is tissue from mouse (or set of pooled mice)
– increase power with more mice on distinct chips
– think of experiment for a single mRNA at a time
© BS Yandell
Experimental Biology 2002
9
Are typical statistical assumptions
warranted for microarrays?
• independence: address at design phase
– want chips independent, but gene spots on chip?
– often expect genes to correlate--coordinated expression
• equal variance
– log (almost) takes care of this--or does it?
– what affects variance? abundance? gene function?
• normality (symmetry, bell-shaped histogram)
– log (almost) transforms to symmetry?
© BS Yandell
Experimental Biology 2002
10
To log or not to log?
• log is natural choice
– tremendous dynamic range (100-1000 fold common)
• intuitive appeal, e.g. concentrations of chemicals (pH)
• fold changes becomes additive
– nice statistical properties ideally
• noise variance roughly constant(?)
• histogram roughly symmetric/normal
– but adjusted values  = PM – MM may be negative
• approximate log transform: normal scores
– there is an exact transform to normality
• close to log() but exact form unknown: -1(F())
• handles negative background-adjusted values
– close approximation easy to compute: X = -1(Fn())
– plot using anti-log to approximate fold changes
© BS Yandell
Experimental Biology 2002
11
Approximate log transform:
normal scores
-4000
0
0.0 0.2 0.4 0.6
Density
0.0008
0.0000
•  = PM – MM
• log10()
note dropped data
proportion
whole chip
23% dropped
proportion
5% dropped
4000
0
PM-MM
1
2
3
4
log10(PM-MM)
0.35
0.05
1 6 12 19 26
rank(PM-MM)
© BS Yandell
0.20
proportion
0.035
0.020
• rank()
• X = -1(Fn())
squish blocks into
bell shaped curve
proportion
sample of 30
Experimental Biology 2002
-2 -1 0 1 2
X = normal score
12
How do we analyze multiple genes?
• assume transformed expression is roughly normal
– at least roughly symmetric
– or use methods that account for data shape
• find common patterns of differential expression
– compare genes across conditions
– how can we combine gene patterns?
• use common patterns in noise
– is variation in noise constant? probably not
• mRNA abundance, gene function, gene-to-gene variability
– how to model changing variation easily?
• let design drive analysis
– linear model based on experimental design
– incorporate sources of variation
© BS Yandell
Experimental Biology 2002
13
Gene-specific model for data analysis
• fit linear model with conditions, genes, replicates
– Xcgr =  + Cc + Gg + Dcg + Ncgr
•
•
•
•
c = condition; g = gene; r = replicate
Cc = 0 if arrays normalized separately
Dcg = differential expression for condition c, gene g
Kerr Churchill (2001)
• mean abundance of gene g: Ag = X•g•
• differential expression:
Dg = D1g – D2g
– contrast among conditions = “signal”
– lean vs. obese, B6 vs. BTBR, …
© BS Yandell
Experimental Biology 2002
14
How to assess differential expression?
• differential expression:
Dg = c wc Xcg•
– Dg =  wc Dcg +  wc Ncg•
– Var(Dg) = g2 + g2/R = signal + noise
• standardized contrasts:  wc = 0,  wc2 = 1
• gene-specific variance of difference
– Var(Dg) = g2/R
– Var(Dg) = g2 + g2/R
no differential expression
differential expression
• infer gene-specific differential expression
– is signal g “large” relative to noise g ?
– how to estimate SDg = g ?
© BS Yandell
Experimental Biology 2002
15
Two ways to measure noise SD
• SD decreases with abundance
– mechanics of hybridization, reading
– SDg2 = (Ag)2
• can estimate without replication
• combine information across genes
– Newton et al. (2001), Roberts et al. (2000), Lin et al. (2001)
• SD varies from gene to gene
– biochemistry of specific mRNA
– SDg2 = gene-specific g2
• need “substantial” replication (say 5?)
• analyze genes separately
– Efron et al. (2002), Lönnstedt Speed (2001)
© BS Yandell
Experimental Biology 2002
16
using replicates
abundance-based SD
using mean contrasts
95% 82 limits
0.05
knockout
SD = spread
0.10 0.20
0.50
gene-specific SD
1.00
2.00
Are SDs constant across genes?
0.02
wild type
0.05
© BS Yandell
0.20
Experimental Biology 2002
0.50
2.00 5.00
average intensity
20.00
17
How to Estimate Spread of Noise?
• focus on genes with no differential expression
– assume SD changes with abundance Ag
– use robust estimate SDg = (Ag) across genes
– screens out changing genes as “outliers”
• focus on replication
– measure expression noise by deviations from mean
• SDg2 = cr(Xcgr – Xcg• )2 / 1
• Combine ideas into gene-specific hybrid
– Gene-specific SDs vary around (Ag)
• “prior” g2 ~ inv-2(0, (Ag)2)
– combines two “statistically independent” estimates
© BS Yandell
Experimental Biology 2002
18
Adipose: no replication
low abundance
5.00
(some  < 0)
obese vs. lean
0.20 0.50
2.00
Clustering
(Nadler et al. 2000)
significant change
B6
F1
BTBR
0.05
lean
obese
0.02
© BS Yandell
0.10
0.50
2.00
average intensity
10.00
Experimental Biology 2002
19
Why Worry about Low Abundance Genes?
• expression may be at or below background level
– background adjustment:  = PM – MM
• removes local “geography”
• allows comparison within and between chips
• can be negative--problem with log transform
– large measurement variability
• early technology (bleeding edge)
• do next generation chips really fix this?
– low abundance genes
• mRNA virtually absent in one condition
• could be important: transcription factors, receptors, regulators
• high prevalence across genes on a chip
– up to 25% per early Affy chips (reduced to 3-5% with www.dChip.org)
– 10-50% across multiple conditions
• low abundance signal may be very noisy
– 50% false positive rate even after adjusting for variance
– may still be worth pursuing: high risk, high research return
© BS Yandell
Experimental Biology 2002
20
Adipose: What was Found?
• transcription factors
– I-kB modulates transcription - inflammatory processes
– RXR nuclear hormone receptor - forms heterodimers
with several nuclear hormone receptors
• regulatory proteins
– protein kinase A
– glycogen synthase kinase-3
• roughly 100 genes
– 90 new since Nadler (2000) PNAS
– but 50% false positives!
B6
F1
BTBR
lean
obese
© BS Yandell
Experimental Biology 2002
21
Robust SD varying with abundance
• median & median absolute deviation (MAD)
– robust to outliers (e.g. changing genes)
– easy to compute
– adapt to patterns in data rather than idealized model
• partition genes into slices based on abundance Ag
– use many slices to assess how SD varies
– ~30 genes per slice for Affy mouse chips (400 slices)
• smooth median & MAD over slices
– automated smoothing splines (Wahba 1990)
– smoothes out slice-to-slice chatter
© BS Yandell
Experimental Biology 2002
22
Median & MAD by abundance slices
(1% sample)
MAD = median absolute deviation
0.0
D = differential expression
-2.0 -1.5 -1.0 -0.5 0.0 0.5
absolute deviation from median
0.5
1.0
1.5
1.0
2.0
median
-1
© BS Yandell
0
1
A = average intensity
2
Experimental Biology 2002
-1
0
1
A = average intensity
2
23
Adipose: center and spread
low abundance
5.00
(some  < 0)
obese vs. lean
0.20 0.50
2.00
Clustering
(Nadler et al. 2000)
significant change
(Bonferroni limits)
B6
F1
BTBR
0.05
lean
obese
0.02
© BS Yandell
0.10
0.50
2.00
average intensity
10.00
Experimental Biology 2002
24
10
Standardized Genotype Effects
standardized effects
Dg / (Ag)
Obesity Dominance
0
5
add = BTBR - B6
dom = F1-(BTBR+B6)/2
significant change
Bonferroni lines
B6
F1
BTBR
-5
lean
obese
-10
© BS Yandell
-5
0
5
Experimental Biology 2002
Obesity Additive
25
Are strain differences real?
SREBP1
BTBR
SCD1
2.6
2.3
B6
BTBR
BTBR
Experimental Biology 2002
BTBR
G6Pase
0.0
1.0
2.0
1.0 1.4 1.8 2.2
2.2
1.8
1.4
B6
B6
3.0
PEPCK
BTBR
© BS Yandell
BTBR
2.9
2
-2 -1 0
BTBR
islet
BTBR
B6
PPARalpha
ACO
B6
muscle
BTBR
1
0.5
B6
B6
B6
PPARgamma
-0.5 0.0
few d.f. per gene
Can we trust SDg ?
1.2
1.9
1.7
-2.5
B6
noise negligible?
liver
GPAT
1.6
-1.5
similar pattern
parallel lines
no interaction
fat
FAS
2.1
strain differences?
B6
BTBR
B6
BTBR
26
Improving on gene-specific SD
• gene-specific SD from replication
– SDg2 = cr(Xcgr – Xcg• )2 / 1
• robust abundanced-based estimate
– (Ag) = smoothed MAD
• Combine ideas into gene-specific hybrid
– “prior” g2 ~ inv-2(0, (Ag)2)
– “posterior” shrinkage estimate
1SDg2 + 0(Ag)2
1 + 0
– combines two “statistically independent” estimates
© BS Yandell
Experimental Biology 2002
27
1.00
SD for strain differences
gene-specific g
0.50
smooth of g
main effects
liver (Ag)
interaction
fat-liver (Ag)
B6
0.05
SD = spread
0.10
0.20
fat (Ag)
B6
fat
liver
muscle
islet
BTBR
BTBR
© BS Yandell
-2
-1
0
1
average intensity
Experimental Biology 2002
2
3
28
B6
0.05
95% 82 limits
new (shrunk) g
size of shrinkage
1g2 + 0(Ag)2
1 + 0
SD = spread
0.10
0.20
gene-specific g
abundance (Ag)
0.50
1.00
Shrinkage Estimates of SD
B6
fat
liver
muscle
islet
BTBR
BTBR
© BS Yandell
-2
Experimental Biology 2002
-1
0
1
average intensity
2
3
29
How good is shrinkage model?
0.8
prior for gene-specific
histogram of ratio
g2 / (Ag)2
empirical Bayes estimates
0 = 5.45,  = 1
2 approximation with 
0 = 3.56,  = .809
used 0 = 8,  = 1 in figures
© BS Yandell
0.0
0.2
2 approximation
0.4
Density
0.6
g2 ~ inv-2(0, (Ag)2)
fudge  to adjust mean
1g2 + 0(Ag)2
1 + 0
0
2
Experimental Biology 2002
4
6
8
10
30
B6
B6
liver
muscle
10
5
0
-5
islet
BTBR
-15
fat
-10
fat-liver interaction
shrinkage-based
abundance-based
9 genes identified
S = (D-center)/spread
15
Effect of SD Shrinkage on Detection
0.02
BTBR
© BS Yandell
0.10
0.50
2.00
10.00
A = average intensity
Experimental Biology 2002
31
-10
-5
significant change
9 genes identified
Bonferroni lines
Fat effect
0
5
10
15
Liver vs. Fat effects
-10
© BS Yandell
-5
Experimental Biology 2002
0
5
Liver effect
10
15
32
How to detect patterns of expression?
• differential expression--or not?
– Dg / (Ag) ~ Normal(0,1) ?
• no differential expression (most genes)
• differential expression more dispersed than N(0,1)
– evaluation of differential expression
• formal test of outliers: multiple comparisons
• posterior probability in differential group
• want to control false positives & false negatives
• general pattern recognition
– in which group does gene belong?
• clustering, discrimination & other multivariate approaches
– linear discriminants are natural extension of ideas here
– are these groups different?
• comparison of functional groups
© BS Yandell
Experimental Biology 2002
33
Multiple Comparisons: a concern?
• many tests performed at once
– goal: detect genes with “large” differential expression
– formality: is Dg / (Ag) ~ Normal(0,1) ?
– practice: use multiple comparisons as guideline
• simple multiple comparisons approach
– Zidak/Bonferroni corrected p-values: p = p1/n
– 13,000 genes with an overall level p = 0.05
• each gene should be tested at level p1 = 1.95*10-6
• differential expression if Dg / (Ag) > 4.62
• is this too conservative? (Dudoit et al. 2000)
– re-envigorated multiple comparisons “industry”
© BS Yandell
Experimental Biology 2002
34
1 e-04
Holms
1 e-08
Sidak
Bonferroni
1 e-12
raw p-value
uniform g/(1+n)
p-value
nominal .05
Holms
Sidak 
Bonferroni
0.05
1 e+00
all multiple comparisons similar
B6
fat
liver
muscle
islet
BTBR
1 e-16
B6
0
2000
4000
8000
10000
12000
rank order
BTBR
© BS Yandell
6000
Experimental Biology 2002
35
pattern of standardized differences
standardized differences
0.2
0.1
probability
Relative Frequency
standard normal
differential expression
Bonferroni cutoff
after Efron et al. (2001)
0.3
Dg / (Ag)
B6
B6
liver
muscle
islet
BTBR
BTBR
© BS Yandell
0.0
fat
-10
-5 -2.9
0
2.8
5
S = (D-center)/spread
Experimental Biology 2002
36
m
4
p
m
7
p
© BS Yandell
10
8
7
6
5
4
2
p
m
5
p
m
8
p
Experimental Biology 2002
m
3
p
m
6
p
m
9
p
9
4
5
6
7
8
9
10
4
5
6
7
8
9
8
4
5
6
7
8
9
10
4
5
6
7
8
7
4
5
6
7
wild type
8
9
10
knockout
4
5
6
relative to gene-to-gene variation
m
10
p
10
1
10
m
9
up or down regulation?
9
10
9
7
5
4
4
5
6
115 significant genes
5-20 genes/group
dropped unknowns
6
7
8
9 functional groups
8
9
10
Comparing gene function groups
37
Related Literature
• comparing two conditions
– log normal: var=c(mean)2
• ratio-based (Chen et al. 1997)
• error model (Roberts et al. 2000; Hughes et al. 2000)
• empirical Bayes (Efron et al. 2002; Lönnstedt Speed 2001)
– gene-specific Dg ~ , var(Dg) ~-1, Zg ~ Bin(p)
– gamma
• Bayes (Newton et al. 2001, Tsodikov et al. 2000)
– gene-specific Xg ~, Zg ~ Bin(p)
• anova (Kerr et al. 2000, Dudoit et al. 2000)
– log normal: var=c(mean)2
– handles multiple conditions in anova model
– SAS implementation (Wolfinger et al. 2001)
• See www.stat.wisc.edu/~yandell/statgen References
© BS Yandell
Experimental Biology 2002
38
R Software Implementation
• quality of scientific collaboration
– hands on experience to researcher
– focus on graphical information content
• needs of implementation
–
–
–
–
quick and visual
easy to use (GUI=Graphical User Interface)
defensible to other scientists
open source in public domain?
• www.r-project.org
– www.bioconductor.org
© BS Yandell
Experimental Biology 2002
39
library(pickgene)
### R library
library(pickgene)
### create differential expression plot(s)
result <- pickgene( data, geneID = probes,
renorm = sqrt(2), rankbased = T )
### print results for significant genes
print( result$pick[[1]] )
### density plot of standardized differences
pickedhist( result, p1 = .05, bw = NULL )
© BS Yandell
Experimental Biology 2002
40
Mapping Gene Expression
as a Quantitative Trait?
• gene expression in segregating population
– assume one gene locus (QTL) influences expression
– create backcross (BC) or intercross (F2)
– map QTL using expression as quantitative trait
• scan entire genome for possible QTL
• MapMaker, QTL Cart or other package
– gene expression may be controlled by other QTL
• multiple genes influenced by same QTL?
– is QTL at a regulatory gene?
• multiple QTL affecting some regulatory gene?
© BS Yandell
Experimental Biology 2002
41
From genes to regulatory genes
• X = expression data from chips for F2 population
– too many gene expressions to map separately
– reduce dimension using multivariate approach
– principle components (singular value decomposition)
• X = UDVT
• V has eigen-genes as rows, individuals as columns
• V = combined expression of coordinated genes
– map first few important rows of V as quantitative traits
– suppose coordination due to gene regulation
• elicit biochemical pathways (Henderson et al. Hoeschele 2001)
– increase power to detect expression-modifying QTL
© BS Yandell
Experimental Biology 2002
42
X = UDVT
© BS Yandell
Alter et al. (2000 PNAS)
yeast cell cycle
Experimental Biology 2002
43
Idea of test for QTL at one locus
(using graphs from West et al. 2000)
colored correlation matrix
genes
X = gene signal
regulators
t-tests at QTL locus
V = regulator signal
© BS Yandell
Experimental Biology 2002
44
Multiple QTLs
• mapping principle component as quantitative trait
– Liu et al. (1996); Zeng et al. (2000)
– multiple interval mapping with interactions
• research groups working on expression QTLs
– Doerge et al. (Purdue)
– Jansen et al. (Waginingen)
• multiple QTL literature
– multiple interval mapping
• Zeng, Kao, et al. (1999, 2000)
– Bayesian interval mapping
• Satagopan et al. (1996); Satagopan, Yandell (1996); Stevens,
Fisch (1998); Silanpää, Arjas (1998, 1999)
© BS Yandell
Experimental Biology 2002
45
Summary
• Why design microarray experiments? (Kerr Churchill)
– chips and samples are expensive: use resources well
– design experiment for one gene with true replication
• Are typical statistical assumptions warranted?
– not automatically--plot your data!
– find transform to symmetry (near normal)
– examine how SD changes with abundance
• How do we combine data analysis across multiple genes?
– keep low abundance data & allow model noise with abundance
– use formal tests as guide to false positive rate
• How can we map gene expression?
– use multivariate summaries to capture functional patterns
– expression may be controlled by other (regulatory) gene
• Ongoing collaboration requires continual dialog
© BS Yandell
Experimental Biology 2002
46
Collaborators
www.stat.wisc.edu/~yandell/statgen
Alan D. Attie3
Hong Lan3
Samuel T. Nadler3
Yi Lin1
Yang Song1
Fei Zou5
Christina Kendziorski4
3UW-Madison
1UW-Madison
© BS Yandell
Biochemistry
Statistics
4UW-Madison Biostatistics
5UNC Biostatistics
Experimental Biology 2002
47