Vraag:
Inzicht in DESeq2-ontwerp, contrast en resultaten
bli
2017-05-19 19:40:14 UTC
view on stackexchange narkive permalink

Ik heb een reeks high-troughput-experimenten met twee genotypen ('WT' en 'prg1') en drie behandelingen ('RT', 'HS30' en 'HS30RT120'), en er zijn twee replica's voor elk van de genotype x behandelingscombinaties.

De gelezen counts voor de genen worden samengevat in een bestand dat ik als volgt laad in R:

  > counts_data <- read.table (" path / to / my / file ", header = TRUE, row.names =" gen ") > colnames (counts_data) [1]" WT_RT_1 "" WT_HS30_1 "" WT_HS30RT120_1 "" prg1_RT_1 "[5]" prg1_HS30_1 "" prg120_HS30_1 " "WT_RT_2" "WT_HS30_2" [9] "WT_HS30RT120_2" "prg1_RT_2" "prg1_HS30_2" "prg1_HS30RT120_2"  

Ik beschrijf de experimenten als volgt:

  > col_data <- DataFrame (geno = c (rep ("WT", times = 3), rep ("prg1", times = 3), rep ("WT", times = 3), rep ("prg1", times = 3)), treat = rep (c ("RT", "HS30", "HS30RT120"), times = 4), rep = c (rep ("1", times = 6), rep ("2", times = 6) ) Row.names = colnames (counts_data)) > col_dataDataFrame met 12 rijen en 3 kolommen geno behandelen rep <character> <character> <character>WT_RT_1 WT RT 1WT_HS30_1 WT WT HS30 1WT_HS30RT120_1 HS30RT120 1prg1_RT_1 PRG1 RT 1prg1_HS30_1 PRG1 HS30 1 ... ... .... ..WT_HS30_2 WT HS30 2WT_HS30RT120_2 WT HS30RT120 2prg1_RT_2 prg1 RT 2prg1_HS30_2 prg1 HS30 2prg1_HS30RT120_2 prg1 HS30RT120 2  

Ik wil een DESeul2-object bouwen dat ik kan gebruiken:

  • vind differentieel tot expressie gebrachte genen wanneer de behandeling varieert voor een bepaald vast genotype
  • of:

    • vind differentieel tot expressie gebrachte genen wanneer het genotype varieert voor een bepaalde vaste behandeling

    In het bioconductor-helpforum denk ik dat ik een enigszins vergelijkbare situatie heb gevonden, en ik las het volgende:

    Probeer een ontwerp van ~ genotype + genotype: condition

    Dan heb je een conditie-effect voor elk genotype-niveau, inclusief het referentieniveau .

    U kunt paren ervan vergelijken met behulp van de lijststijl van het 'contrast'-argument.

    Dit verklaart echter niet hoe u deze' lijststijl 'toepast "naar het" contrast "-argument. En de bovenstaande situatie lijkt asymmetrisch te zijn. Daarmee bedoel ik dat genotype en conditie geen onderling verwisselbare rol lijken te hebben.

    Dus probeerde ik de volgende meer symmetrische formule:

      > dds <- DESeqDataSetFromMatrix (countData = counts_data, colData = col_data, design = ~ geno + treat + geno: treat) > dds <- DESeq (dds)  

    Kan ik nu bijvoorbeeld de differentiële uitdrukking krijgen resultaten bij vergelijking van behandeling "HS30" met "RT" als referentie, in genotype "prg1" ?

    En hoe?

    Als ik het goed begrijp, is het bovenstaande genoemde "lijststijl" gebruikt namen gegeven door de functie resultsNames . In mijn geval heb ik het volgende:

      > resultsNames (dds) [1] "Intercept" "geno_WT_vs_prg1" [3] "treat_HS30RT120_vs_HS30" "treat_RT_vs_HS30" [5] "genoWT.treatHS30RT120" "genoWT.treatRT"  

    Ik denk dat ik een contrast nodig zou hebben tussen "genoprg1.treatRT" en een "genoprg1.treatHS30", maar deze staan ​​niet in de bovenstaande resultatennamen.

    Ik ben verdwaald.

    dit is een meer vergelijkbaar voorbeeld (denk ik): https://support.bioconductor.org/p/63201/. De conclusie is vergelijkbaar met het antwoord van Devon - gebruik LRT.
    Twee antwoorden:
    #1
    +13
    Devon Ryan
    2017-05-19 19:51:01 UTC
    view on stackexchange narkive permalink

    De eenvoudigste manier is om geen wald-test te gebruiken, maar eerder een LRT met een gereduceerd model zonder de factor van belang:

      dds = DESeq (dds, test = "LRT" gereduceerd = ~ geno + geno: behandeling)  

    Het bovenstaande geeft u resultaten voor behandeling ongeacht het niveau, terwijl u nog steeds rekening houdt met een mogelijke interactie (dwz een 'hoofdeffect van de behandeling, ongeacht de type behandeling ').

    Even terzijde, dit is waarschijnlijk een geval waarin de edgeR-voorkeurmanier om groepen genotype-behandelingscombinaties te maken en vervolgens een model van ~ 0 + groep kan uw leven een beetje gemakkelijker maken. U krijgt hoe dan ook dezelfde resultaten (min of meer), maar het zal waarschijnlijk gemakkelijker voor u zijn om in die termen te denken dan te onthouden dat het basisniveau behandeling HS30 en zal zijn geno prg1 .

    Ik heb het volgende geprobeerd: `dds <- DESeqDataSetFromMatrix (countData = counts_data, colData = col_data, design = ~ geno + treat + geno: treat)` dan `dds = DESeq (dds, test =" LRT ", gereduceerd = ~ geno + geno: treat) `maar dit mislukt met` Fout in nbinomLRT (object, full = volledig, gereduceerd = gereduceerd, betaPrior = betaPrior,: minder dan één vrijheidsgraad, misschien staan ​​volledige en gereduceerde modellen niet in de juiste volgorde '. , wat ik wil is niet "ongeacht het type behandeling". Ik wil resultaten voor een bepaalde vaste behandeling of voor een bepaald vast genotype. Misschien was mijn vraag niet duidelijk.
    Ik zou de modelmatrix moeten doornemen om te zien waarom deze de fout veroorzaakt. Hoe dan ook, ik zou willen voorstellen om de route te volgen die de edgeR-auteurs verkiezen en de combinaties van behandeling en genotype in individuele groepen te maken. Dan zijn de contrasten duidelijker.
    "De eenvoudigste manier is om geen wald-test te gebruiken, maar eerder een LRT" - zou je kunnen uitleggen waarom?
    Omdat de behandeling meerdere niveaus heeft, zou men met een wald-test de contrasten op de juiste manier moeten behandelen, terwijl men met een LRT tegelijkertijd alle niveaus in één keer kan weglaten.
    #2
    +8
    bli
    2017-05-22 17:05:18 UTC
    view on stackexchange narkive permalink

    Het lijkt erop dat de truc van "factoren combineren" beschreven in deel 3.3 van het huidige "vignet" van DESeq2 (vanaf mei 2017) onder de titel "Interactie" een manier is om toegang te krijgen tot de gewenste contrasten.

    Het lijkt mogelijk om het direct te doen bij het bouwen van de colData en bij het aanroepen van DESeqDataSetFromMatrix:

    Laten we een gecombineerde "geno" en "treat" factoren toevoegen naar de toekomstige parameter colData :

      > col_data $ geno_treat <- as.factor (plak (col_data $ geno, col_data $ treat, sep = "_")) > col_dataDataFrame 12 rijen en 4 kolommen geno behandelen rep geno_treat <character> <character> <character> <factor>WT_RT_1 WT RT1 WT_RTWT_HS30_1 WT HS30 1 WT_HS30WT_HS30RT120_1 WT HS30RT120 1 WT_HS30RT120prg1_RT_1 PRG1 RT1 prg1_RTprg1_HS30_1 PRG1 HS30 1 prg1_HS30 ... ... ... ... ... WT_HS30_2 WT HS30 2 WT_HS30WT_HS30RT120_2 WT HS30RT120 2 WT_HS30RT120prg1_RT_2 PRG1 RT 2 prg1_RTprg1_HS30_2 PRG1 HS30 2 prg1_HS30prg1_HS30RT120_2 PRG1 HS30RT120 2 prg1_HS30RT120  

    We kunnen nu een ontwerp gebruiken waarbij differentiële expressie wordt verklaard door deze gecombineerde factoren:

      > dds <- DESeqDataSetFromMatrix (countData = counts_data, colData = col_data, design = ~ geno_treat)  

    We voeren de analyse uit:

      > dds <- DESeq (dds)  

    Vervolgens kunnen we resultaten opvragen voor een bepaald contrast tussen dergelijke factorcombinaties. Om bijvoorbeeld de resultaten te hebben voor het effect van behandeling "HS30" tegen de referentietoestand "RT" in genotype "prg1":

      res <- resultaten (dds, contrast = c (" geno_treat "," prg1_HS30 "," prg1_RT "))  


    Deze Q&A is automatisch vertaald vanuit de Engelse taal.De originele inhoud is beschikbaar op stackexchange, waarvoor we bedanken voor de cc by-sa 3.0-licentie waaronder het wordt gedistribueerd.
    Loading...