Statistische Analyse von hochdimensionalen Daten in der ......Einführung RNA Expression...

Transcript of Statistische Analyse von hochdimensionalen Daten in der ......Einführung RNA Expression...

Einführung RNA Expression Statistisches Testen Multiples Testen GWAS

Statistische Analyse vonhochdimensionalen Daten

in der Bioinformatik

Florian Frommlet

Institut für medizinische Statistik,Medizinische Universität Wien

Wien, November 2013

Einführung

DNA Molekül

Zwei komplementäre Stränge bildenDoppelhelix

Vier Basen

• Adenin• Thymin• Cytosin• Guanin

http://neutronsforbreakfast.wordpress.com/

Einführung

DNA Molekül

Zwei komplementäre Stränge bildenDoppelhelix

Vier Basen

• Adenin• Thymin• Cytosin• Guanin

Genetische Information als String inAlphabet mit vier Buchstaben http://neutronsforbreakfast.wordpress.com/

gaacgaatca ttgcaaagag ccaaagatcc aaaatttgca acaaaaacaa aaactctacc

Welche genetische Information?

Codierung von Proteinen

• Makromoleküle aus 20Aminosäuren

• Grundbausteine aller Zellen• Unglaubliche Vielfalt anAufgaben

Codierung von RNARibonukleinsäure - Umsetzung vongenetischer Information, aber auchviele andere Aufgaben

http://techglimpse.com

Welche genetische Information?

Codierung von Proteinen

• Makromoleküle aus 20Aminosäuren

• Grundbausteine aller Zellen• Unglaubliche Vielfalt anAufgaben

Codierung von RNARibonukleinsäure - Umsetzung vongenetischer Information, aber auchviele andere Aufgaben

http://techglimpse.com

Was ist ein Gen?

Übliche Definition:Abschnitt der DNA, der ein ProteincodiertTriplets von DNA - Basencodieren Aminosäuren⇒ 64 Möglichkeiten (Redundanz)

http://de.academic.ru

Was ist ein Gen?

Übliche Definition:Abschnitt der DNA, der ein ProteincodiertTriplets von DNA - Basencodieren Aminosäuren⇒ 64 Möglichkeiten (Redundanz)

http://de.academic.ru

Standardmodell der Genetik

http://kvhs.nbed.nb.ca

Technologien zum Messen von RNA - Expression

Älter: Microarrays

• Chip mit kurzen DNA Stückenvon Genen an Spots(jeweils nur ein Strang)

• Hybridisierung: Anlagerungvon komplementärem DNAoder RNA Strang

Neuer: RNA-SeqKomplexere Technologie, die erlaubtmittels next generation sequencingeinen generellen Überblick über diein einer Zelle vorhandene RNA zuerhalten http://pragmatyczny17.siam.im

Technologien zum Messen von RNA - Expression

Älter: Microarrays

• Chip mit kurzen DNA Stückenvon Genen an Spots(jeweils nur ein Strang)

• Hybridisierung: Anlagerungvon komplementärem DNAoder RNA Strang

Neuer: RNA-SeqKomplexere Technologie, die erlaubtmittels next generation sequencingeinen generellen Überblick über diein einer Zelle vorhandene RNA zuerhalten http://pragmatyczny17.siam.im

Prinzip von Micro Arrays

Was wird gemessen?DNA Proben mit FluorophorengekennzeichnetJe mehr Hybridisierung an einemSpot desto stärker das Farbsignal

Zweifarben MicroarraysZwei Gruppen mit unterschiedlichenFarben gekennzeichnet (z. Bsp. rotund grün)Gelb: Beide Gruppen exprimiertDunkel: Keine von beiden

http://pragmatyczny17.siam.im

Prinzip von Micro Arrays

Was wird gemessen?DNA Proben mit FluorophorengekennzeichnetJe mehr Hybridisierung an einemSpot desto stärker das Farbsignal

Zweifarben MicroarraysZwei Gruppen mit unterschiedlichenFarben gekennzeichnet (z. Bsp. rotund grün)Gelb: Beide Gruppen exprimiertDunkel: Keine von beiden

http://pragmatyczny17.siam.im

Datenaufbereitung

Vor der eigentlichen Datenanalyse müssen die Bilddaten aufgearbeitetwerden ⇒ Eigene Vorlesung

Wesentliche Schritte:• Bildanalyse (speziell interessant bei Zweifarben-Arrays)• Qualitätskontrolle (Viele Fehlermöglichkeiten)• Transformation und Normalisierung (um danach einfache statistischeVerfahren anwenden zu können)

• Behandlung von fehlenden Werten• etc.

Am Ende dieser Prozedur erhält man bei Zweifarben-Arrays für jedenSpot eine Zahl, die mißt ob es für das jeweilige Gen einen Unterschied inder Genexpression zwischen den beiden jeweiligen Proben gibt.

Datenaufbereitung





Datenaufbereitung





Datenaufbereitung





Datenaufbereitung





Datenaufbereitung





Beispiel: Studie über Prostatakrebs

Singh et al. (2002)“Gene expression correlates ofclinical prostate cancer behavior,”Cancer Cell 1, pp. 203-209.

Microarray Experiment: Von jedemPatien gesunde und Krebszellen

• Gleason Score als Maß fürDifferenzierungsgrad der Zelle

• 52 Individuen: 26 davon mitschlechter Differenzierung

• etwa 12000 GeneWikipedia

Frage: Gibt es Gene die Differenzierungsgrad beeinflussen?

Studie über ProstatakrebsEin spezielles GenHomo sapiens mRNA for RET finger protein-like 3

Graphische Darstellung: Boxplot

Vergleich für ein einziges Gen

Nullhypothese: Expression in beiden Gruppen ist gleich

H0 : µ1 = µ2

Alternative: Expression in beiden Gruppen ist verschieden

HA : µ1 6= µ2

Vergleiche Mittelwerte zwischen den beiden Gruppen

Gruppe 1: x̄1 = 17.0769, Gruppe 2: x̄2 = 8.8846

Je mehr die Daten streuen desto weniger Aussagekräftig ist dieserUnterschied zwischen den beiden Gruppen



H0 : µ1 = µ2


HA : µ1 6= µ2


Gruppe 1: x̄1 = 17.0769, Gruppe 2: x̄2 = 8.8846




H0 : µ1 = µ2


HA : µ1 6= µ2


Gruppe 1: x̄1 = 17.0769, Gruppe 2: x̄2 = 8.8846


Der Zweistichproben t-Test

t-Test Statistik T :

T =x̄1 − x̄2

S

wobei S2 geeigneter Schätzer für die Varianz der Mittelwertsdifferenz.

TestentscheidungFalls |T | größer als kritischer Wert ⇒ Entscheidung für HASonst Beibehaltung von H0

Der Zweistichproben t-Test

t-Test Statistik T :

T =x̄1 − x̄2

S

wobei S2 geeigneter Schätzer für die Varianz der Mittelwertsdifferenz.

TestentscheidungFalls |T | größer als kritischer Wert ⇒ Entscheidung für HASonst Beibehaltung von H0

Wahrscheinlichkeitsverteilung von TTheoretische Verteilung von T = x̄1−x̄2S unter Nullhypothese:t-verteilt mit n − 2 Freiheitsgraden

Dichte der t-Verteilung mit df Freiheitsgraden

Bestimmung des kritischen WertsUnter Annahme dass H0 stimmt, suche symmetrischen Bereich wo H0mit Wahrscheinlichkeit α abgelehnt wirdd.h. Nullhypothese wird fälschlich verworfen (Fehler 1. Art)

Wahl des kritischen Werts für α = 0.05 und df = 50

t-Test für unser spezielles Gen

x1 = 17.08, x2 = 8.88S = 1.9249⇒ T = 4.26

Kritischer Wert: 2.0086

⇒ Entscheidung für H1

p-Wert: Wahrscheinlichkeit unter H0 den Wert T oder noch einenextremeren Wert zu beobachten

|T | > kritischer Wert ⇔ p-Wert < α

Im Beispiel: p = 9.1372 ∗ 10−5


x1 = 17.08, x2 = 8.88S = 1.9249⇒ T = 4.26





Im Beispiel: p = 9.1372 ∗ 10−5


x1 = 17.08, x2 = 8.88S = 1.9249⇒ T = 4.26





Im Beispiel: p = 9.1372 ∗ 10−5

Prinzip des statistischen Testens

Entscheidung z. Bsp. basierend auf t-Test Statistik T :

|T | > kritischer Wert ⇒ H0 abgelehnt|T | ≤ kritischer Wert ⇒ H0 beibehalten

H0 beibehalten H0 abgelehntH0 wahr Fehler 1. ArtH1 wahr Fehler 2. Art

Kontrolliere Fehler 1. Art zum Niveau α (z.B. α = 0.05, oder α = 0.01)

PowerWahrscheinlichkeit einen tatsächlichen Unterschied zu erkennen

Je kleiner α ⇒ desto kleiner die Power




H0 beibehalten H0 abgelehntH0 wahr Fehler 1. ArtH1 wahr Fehler 2. Art Power




Multiples Testproblem

Ein statistischer Test: Wahrscheinlichkeit für falsche Entdeckung ist α

Was passiert bei vielen (m) Tests?Family wise error rate:

FWER := P(Mindestens eine falsche Entdeckung)

In unserem Beispiel: m = 12000 GeneWie groß ist FWER für α = 0.05 ?





In unserem Beispiel: m = 12000 GeneWie groß ist FWER für α = 0.05 ?Unter der (gewagten) Annahme, dass die Tests unabhängig sind:

FWER = 1− (1− α)m = 1− 4.8224 ∗ 10−268

FWER praktisch nicht von 1 unterscheidbar

Man erwartet 12000 ∗ 0.05 = 600 falsche Entdeckungen


Es bedarf einer Korrektur für multiples Testen

Bonferroni Korrektur

Klassische (1936) und extrem einfache Korrekturmethode:Teste zum Signifikanzniveau αBon = α/mBeispiel:

α = 0.05,m = 10 ⇒ αBon = 0.005α = 0.05,m = 100 ⇒ αBon = 0.0005

Bonferroni kontrolliert FWER weil

P

m⋃j=1

Ai

≤ m∑j=1

P(Ai )

Bonferroni Korrektur sehr populär, aber Problem mit geringer Power

In unserem Beispiel αBon = 4.1667 ∗ 10−6⇒ Gen mit p-Wert 9.1372 ∗ 10−5 nicht mehr signifikant



α = 0.05,m = 10 ⇒ αBon = 0.005α = 0.05,m = 100 ⇒ αBon = 0.0005


P

m⋃j=1

Ai

≤ m∑j=1

P(Ai )





α = 0.05,m = 10 ⇒ αBon = 0.005α = 0.05,m = 100 ⇒ αBon = 0.0005


P

m⋃j=1

Ai

≤ m∑j=1

P(Ai )



False Discovery Rate

Konzept der FDREingeführt von Benjamini und Hochberg (1995)

FDR = E

(V

R

)R . . . Gesamtanzahl der abgelehnten HypothesenV . . . Anzahl der falschen EntdeckungenV /R = 0 falls R = 0

FDR: Erwarteter relativer Anteil an falschen Entdeckungen

FWER zu kontrollieren ist eine stärkere Anforderung als FDR,

⇒ größere Power, dafür mehr falsche Entdeckungen



FDR = E

(V

R







FDR = E

(V

R





Benjamini - Hochberg Prozedur

Vorgangsweise

1. Ordne alle p-Werte:p[1] ≤ · · · ≤ p[m]

2. Bestimmek = argmaxi

{p[i ] ≤ iαm

}3. Ablehnung aller Hypothesen mit

p-Wert ≤ p[k]

Beispiel: m = 5p[1] = 0.007 ≤ 0.01p[2] = 0.013 ≤ 0.02p[3] = 0.031 > 0.03p[4] = 0.039 ≤ 0.04p[5] = 0.231 > 0.05

Verwerfe 4 Hypothesen

BH kontrolliert FDR zum Level α


Vorgangsweise



{p[i ] ≤ iαm


p-Wert ≤ p[k]

Beispiel: m = 5p[1] = 0.007 ≤ 0.01p[2] = 0.013 ≤ 0.02p[3] = 0.031 > 0.03p[4] = 0.039 ≤ 0.04p[5] = 0.231 > 0.05




Vorgangsweise



{p[i ] ≤ iαm


p-Wert ≤ p[k]

Beispiel: m = 5p[1] = 0.007 ≤ 0.01p[2] = 0.013 ≤ 0.02p[3] = 0.031 > 0.03p[4] = 0.039 ≤ 0.04p[5] = 0.231 > 0.05




Vorgangsweise



{p[i ] ≤ iαm


p-Wert ≤ p[k]

Beispiel: m = 5p[1] = 0.007 ≤ 0.01p[2] = 0.013 ≤ 0.02p[3] = 0.031 > 0.03p[4] = 0.039 ≤ 0.04p[5] = 0.231 > 0.05




Vorgangsweise



{p[i ] ≤ iαm


p-Wert ≤ p[k]

Beispiel: m = 5p[1] = 0.007 ≤ 0.01p[2] = 0.013 ≤ 0.02p[3] = 0.031 > 0.03p[4] = 0.039 ≤ 0.04p[5] = 0.231 > 0.05



Studie über Prostatakrebs

Anzahl signifikanter Gene

α = 0.05m = 12000

• Keine Korrektur: 1967• Bonferroni: 1• Benjamini Hochberg: 5

10 kleinsten p-Wertep-Wert iαn

p[1] = 3.0 ∗ 10−6 4.17 ∗ 10−6p[2] = 9.5 ∗ 10−6 8.33 ∗ 10−6p[3] = 10.5 ∗ 10−6 12.50 ∗ 10−6p[4] = 11.1 ∗ 10−6 16.67 ∗ 10−6p[5] = 15.2 ∗ 10−6 20.83 ∗ 10−6p[6] = 70.5 ∗ 10−6 25.00 ∗ 10−6p[7] = 91.4 ∗ 10−6 29.17 ∗ 10−6p[8] = 207.8 ∗ 10−6 33.33 ∗ 10−6p[9] = 236.2 ∗ 10−6 37.50 ∗ 10−6p[10] = 296.5 ∗ 10−6 41.67 ∗ 10−6



α = 0.05m = 12000



p[1] = 3.0 ∗ 10−6 4.17 ∗ 10−6p[2] = 9.5 ∗ 10−6 8.33 ∗ 10−6p[3] = 10.5 ∗ 10−6 12.50 ∗ 10−6p[4] = 11.1 ∗ 10−6 16.67 ∗ 10−6p[5] = 15.2 ∗ 10−6 20.83 ∗ 10−6p[6] = 70.5 ∗ 10−6 25.00 ∗ 10−6p[7] = 91.4 ∗ 10−6 29.17 ∗ 10−6p[8] = 207.8 ∗ 10−6 33.33 ∗ 10−6p[9] = 236.2 ∗ 10−6 37.50 ∗ 10−6p[10] = 296.5 ∗ 10−6 41.67 ∗ 10−6



α = 0.05m = 12000



p[1] = 3.0 ∗ 10−6 4.17 ∗ 10−6p[2] = 9.5 ∗ 10−6 8.33 ∗ 10−6p[3] = 10.5 ∗ 10−6 12.50 ∗ 10−6p[4] = 11.1 ∗ 10−6 16.67 ∗ 10−6p[5] = 15.2 ∗ 10−6 20.83 ∗ 10−6p[6] = 70.5 ∗ 10−6 25.00 ∗ 10−6p[7] = 91.4 ∗ 10−6 29.17 ∗ 10−6p[8] = 207.8 ∗ 10−6 33.33 ∗ 10−6p[9] = 236.2 ∗ 10−6 37.50 ∗ 10−6p[10] = 296.5 ∗ 10−6 41.67 ∗ 10−6

Weitere Methoden

Multiples Testen ein sehr aktiver Forschungsbereich

Einige weitere Stichworte• Permutationstests• Baysianische Modellselektion• Empirical Bayes

Literatur für Microarrayanalyse

• Speed (2003) Statistical analysis of gene expression microarray data• Lee (2004) Analysis of microarray gene expression data• Zhang (2006) Advanced analysis of gene expression microarray data• Mallick, Gold (2009) Bayesian analysis of microarray gene expression

data• etc.

Weitere Methoden

Multiples Testen ein sehr aktiver Forschungsbereich

Einige weitere Stichworte• Permutationstests• Baysianische Modellselektion• Empirical Bayes

Literatur für Microarrayanalyse

• Speed (2003) Statistical analysis of gene expression microarray data• Lee (2004) Analysis of microarray gene expression data• Zhang (2006) Advanced analysis of gene expression microarray data• Mallick, Gold (2009) Bayesian analysis of microarray gene expression

data• etc.

Genome Wide Associations Studies

Assoziationsstudie

Suche nach Regionen derDNA die im Zusammenhangmit Merkmalen stehen

Merkmale:• Quantitativ (Größe)• Dichotom (Krankheit)• etc.

Genetischer Marker:Position der DNA wo esUnterschiede zwischenIndividuen gibt http://kvhs.nbed.nb.ca


Assoziationsstudie





Assoziationsstudie




SNPs als genetische Marker

Single NucleotidePolymorphism

SNP: Punktmutation

Beim Menschen fast 20 MillionenSNPs bekannt

HapMap ProjektKarte von SNPs für 270 Individuen

Wikipedia

SNPs als genetische Marker

Single NucleotidePolymorphism

SNP: Punktmutation

Beim Menschen fast 20 MillionenSNPs bekannt

HapMap ProjektKarte von SNPs für 270 Individuen

Wikipedia

SNP ArraysBis zu 1 Million SNPs auf einem Array

Technologie ähnlich zu Microarrays• Zwei Spots für jeweiligeVariante von einem SNP

• Hybridisierung• Markierung mit Fluorophoren

ChromosomenpaarZwei Allele3 mögliche Genotypen: aa, aA, AAHomozygot vs. Heterozygot

Clusteralgorithmen zur Bestimmungdes Genotyps

from BLOG of Bryce Christensen

GWAS Datenstruktur

Y ← X1, . . . ,Xm

• n Individuen typisch n > 1000• m SNPs typisch m > 100000• Y . . . n - Vektor der Merkmale• Xi . . . n - Vektor der GenotypenCodierung z. Bsp Xi ∈ {−1, 0, 1}

Frage:Welche Xi sind mit dem Merkmal assoziiert?

Oft hat SNP nicht unmittelbar selbst Einfluss auf Merkmal sondern istnur Indikator für Gen in der Nähe (Linkage disequilibrium)

GWAS Datenstruktur

Y ← X1, . . . ,Xm

• n Individuen typisch n > 1000• m SNPs typisch m > 100000• Y . . . n - Vektor der Merkmale• Xi . . . n - Vektor der GenotypenCodierung z. Bsp Xi ∈ {−1, 0, 1}

Frage:Welche Xi sind mit dem Merkmal assoziiert?

Oft hat SNP nicht unmittelbar selbst Einfluss auf Merkmal sondern istnur Indikator für Gen in der Nähe (Linkage disequilibrium)

Test von einzelnen Markern

Die einfachste und derzeit auch am weitesten verbreitete Methode zurDatenanalyse von GWAS besteht darin, jeden SNP einzeln zu untersuchen

Test von individuellen MarkernVerschiedene Tests für Zusammenhang zwischen Y und Xi

• Y quantitativ: Lineare Regression, ANOVA, etc.• Y dichotom: χ2-Test, Cochran-Armitage, etc.

In jedem Fall wiederum multiple Testkorrektur notwendig

Test von einzelnen Markern

Die einfachste und derzeit auch am weitesten verbreitete Methode zurDatenanalyse von GWAS besteht darin, jeden SNP einzeln zu untersuchen

Test von individuellen MarkernVerschiedene Tests für Zusammenhang zwischen Y und Xi

• Y quantitativ: Lineare Regression, ANOVA, etc.• Y dichotom: χ2-Test, Cochran-Armitage, etc.

In jedem Fall wiederum multiple Testkorrektur notwendig

Alternative: Modellselektion

Im Falle von quantitativen Merkmalen mittels Regressionsmodell

i. e. Y = β0 + βi1Xi1 + βi2Xi2 + · · ·+ βikXik + �,

Modell mit k SNPs

Es gibt 2m − 1 mögliche ModelleWie findet man das ’richtige’?

Derzeit mein Forschungsbereich (WWTF Projekt)• Modellselektionskriterien für m > n• Suchstrategien in dem riesigen Raum der Modelle• Sparsity hilfreich (k � m) - Entwicklung von Theorie• Effiziente Implementation• etc.




Modell mit k SNPs






Modell mit k SNPs






Modell mit k SNPs






Modell mit k SNPs






Modell mit k SNPs



EinführungRNA ExpressionStatistisches TestenMultiples TestenGWAS

Statistische Analyse von hochdimensionalen Daten in der ......Einführung RNA Expression...

Documents

Transcript of Statistische Analyse von hochdimensionalen Daten in der ......Einführung RNA Expression...