[BioC] Need help: no MTC possible

James W. MacDonald jmacdon at uw.edu
Tue Oct 16 15:48:13 CEST 2012


Hi Suparna,

On 10/16/2012 5:58 AM, suparna mitra wrote:
> Hello group,
> Related to my previous post, I further tried arrayweight as:
>
> > f.invivo <- factor(InVivoTargets$Treatment, levels = c("A", "R", "T"))
>
> > design.invivo <- model.matrix(~0 + f.invivo)
>
> > colnames(design.invivo) <- c("A", "R", "T")
>
> > design.invivo
>
>    A R T
>
> 1  1 0 0
>
> 2  1 0 0
>
> 3  1 0 0
>
> 4  1 0 0
>
> 5  1 0 0
>
> 6  1 0 0
>
> 7  0 1 0
>
> 8  0 1 0
>
> 9  0 1 0
>
> 10 0 1 0
>
> 11 0 1 0
>
> 12 0 1 0
>
> 13 0 0 1
>
> 14 0 0 1
>
> 15 0 0 1
>
> 16 0 0 1
>
> 17 0 0 1
>
> 18 0 0 1
>
> attr(,"assign")
>
> [1] 1 1 1
>
> attr(,"contrasts")
>
> attr(,"contrasts")$f.invivo
>
> [1] "contr.treatment"
>
>
> >
>
> > arrayw <- arrayWeightsSimple(rmaOligoinvivo, design.invivo)
>
> > fit <- lmFit(rmaOligoinvivo, design.invivo, weights=arrayw)
>
> > arrayw
>
>         1         2         3         4         5         6         7 
>         8         9        10        11        12        13        14
>
> 0.3749711 0.8578285 1.9289731 1.2390065 0.8116796 1.7846502 1.0741852 
> 1.4277605 0.6533368 0.7637412 1.2647738 1.4520790 0.8309346 0.9328655
>
>        15        16        17        18
>
> 1.1926458 0.7280477 0.5130294 1.8503073
>
> > contrast.matrix.invivo <- makeContrasts(R-A, T-R, T-A,levels = 
> design.invivo)
>
> > fit2<-contrasts.fit(fit, contrast.matrix.invivo)
>
> > fit2 <- eBayes(fit2)
>

Looks good to me.

Best,

Jim


> >
>
> Can anybody please suggest if I am doing it right? Actually being new 
> in this I am bit afraid to make errors.
> Thanks,
> Suparna.
>
> On 16 October 2012 10:36, suparna mitra <smitra at liverpool.ac.uk 
> <mailto:smitra at liverpool.ac.uk>> wrote:
>
>     Dear James,
>       Thanks for your suggestion. I was reading arrayWeights package
>     in limma.
>     But being novice in bioC I have one confusion. Should I
>     perform arrayWeights on normalized (rmaOligo) expression data or
>     on the raw data?
>
>     This is what i have done so far:
>
>     > sessionInfo()
>     R version 2.15.1 (2012-06-22)
>     Platform: i386-apple-darwin9.8.0/i386 (32-bit)
>
>     locale:
>     [1] en_GB.UTF-8/en_GB.UTF-8/en_GB.UTF-8/C/en_GB.UTF-8/en_GB.UTF-8
>
>     attached base packages:
>     [1] stats     graphics  grDevices utils     datasets  methods   base
>
>     other attached packages:
>      [1] statmod_1.4.15              limma_3.12.1              
>      annotate_1.34.1             hugene10stprobeset.db_8.0.1
>     org.Hs.eg.db_2.7.1
>      [6] BiocInstaller_1.4.7         affycoretools_1.28.0      
>      KEGG.db_2.7.1               GO.db_2.7.1                
>     AnnotationDbi_1.18.1
>     [11] affy_1.34.0                 Biobase_2.16.0            
>      BiocGenerics_0.2.0          pd.hugene.1.0.st.v1_3.6.0  
>     RSQLite_0.11.1
>     [16] DBI_0.2-5                   oligo_1.20.4              
>      oligoClasses_1.18.0
>
>
>      rmaOligoinvivo = oligo::rma(InVivodat1)
>     Background correcting
>     Normalizing
>     Calculating Expression
>
>     > maplot(rmaOligoinvivo)
>     >hist(rmaOligoinvivo)
>     > InVivoTargets=readTargets("~/Desktop/Recent/Liverpool-work-related/Micro_RawData/InVivoTargets.txt")
>     > InVivoTargets
>      FileName Treatment
>     1   MC1       A
>     2   MC2       A
>     3   MC3       A
>     4   MC4       A
>     5   MC5       A
>     6   MC6       A
>     7   MC7       R
>     8   MC8        R
>     9   MC9        R
>     10 MC10        R
>     11 MC11        R
>     12 MC12        R
>     13 MC13       T
>     14 MC14        T
>     15 MC15        T
>     16 MC16        T
>     17 MC17        T
>     18 MC18        T
>
>     f.invivo <- factor(InVivoTargets$Treatment, levels = c("A", "R", "T"))
>
>     design.invivo <- model.matrix(~0 + f.invivo)
>
>     > colnames(design.invivo) <- c("A", "R", "T")
>
>     > fit.invivo <- lmFit(rmaOligoinvivo, design.invivo)
>
>     > contrast.matrix.invivo <- makeContrasts(R-A, T-R, T-A,levels =
>     design.invivo)
>
>     > fit2.invivo <- contrasts.fit(fit.invivo, contrast.matrix.invivo)
>
>     > fit2.invivo <-eBayes(fit2.invivo)
>
>     Thanks a lot,
>     Suparna.
>
>
>     On 15 October 2012 14:33, James W. MacDonald <jmacdon at uw.edu
>     <mailto:jmacdon at uw.edu>> wrote:
>
>         Hi Suparna,
>
>
>         On 10/15/2012 7:01 AM, suparna mitra wrote:
>
>             Hi all,
>                I have been working in a project where I have
>             Affymetrix Hgene 1.0 St V1
>             data. And I have tree groups of patients having 6 samples
>             each. I tried to
>             perform rma normalization and to filter my data based on
>             expression values
>             20%. After that went for unpaired t-test to test each two
>             combination of
>             groups. But the problem is my data is extremely variable.
>             I have tried to filter my genes based on variance and/or
>             CV before testing,
>             to try to reduce the number of genes entering your test
>             and multiple
>             correction.  But with different reasonable filtering also
>             I am with no
>             luck. And I don't have the option to increase sample size
>             of my project.
>             Further I tried to check for the bad samples and bad
>             probes from
>             experimentand remove outlier if these are not of interest.
>             Still the same
>             when run t-test (and other possible test like
>             Mann-Whitney) with MTC there
>             are no genes.
>             On the other hand if I go on with out MTC and select a
>             good p value cutoff
>             and reasonable fold change I get a list of significant
>             gene which may be
>             good or reasonable for my study. but the problem is I
>             somehow need to
>             justify the method for my finding. Do you know any study
>             or paper where
>             anybody has treated their data without MTC?
>             My main concern is if I find a good story matching
>             biological prospective,
>             would it be anyhow possible to justify the method without MTC?
>
>
>         It's not clear to me what you are doing here - when you filter
>         on variance are you keeping or removing the high variability
>         genes (keeping, I hope)? I am also not sure what MTC stands
>         for - is this multiple test correction?
>
>         Anyway, assuming I have things correct, some suggestions.
>         First, you might want to use array weights when fitting your
>         model. If you have a lot of intra-group variability, this will
>         tend to help.
>
>         Second, the t-statistic is the universally most powerful test
>         (assuming the underlying data are relatively hump-shaped), so
>         going to a non-parametric test will usually reduce rather than
>         increase power to detect differences.
>
>         Third, univariate tests are arguably not the most
>         sophisticated way of analyzing expression data, and you might
>         get better (or at least more satisfactory) results if you
>         instead looked at analyzing for groups of genes rather than
>         individually.
>
>         Depending on your experiment, you could accomplish this task
>         with a gene set analysis (there are multiple ways of doing
>         this - perhaps the easiest being romer() and roast() in
>         limma), or if you have phenotypic data, especially continuous
>         measures, a WGCNA analysis might be of some use.
>
>         Best,
>
>         Jim
>
>
>             Thanks a lot,
>             Suparna.
>
>                     [[alternative HTML version deleted]]
>
>             _______________________________________________
>             Bioconductor mailing list
>             Bioconductor at r-project.org <mailto:Bioconductor at r-project.org>
>             https://stat.ethz.ch/mailman/listinfo/bioconductor
>             Search the archives:
>             http://news.gmane.org/gmane.science.biology.informatics.conductor
>
>
>         -- 
>         James W. MacDonald, M.S.
>         Biostatistician
>         University of Washington
>         Environmental and Occupational Health Sciences
>         4225 Roosevelt Way NE, # 100
>         Seattle WA 98105-6099
>
>
>
>
>     -- 
>     Dr. Suparna Mitra
>     Wolfson Centre for Personalised Medicine
>     Department of Molecular and Clinical Pharmacology
>     Institute of Translational Medicine University of Liverpool
>     Block A: Waterhouse Buildings,  L69 3GL Liverpool
>
>     Tel. +44 (0)151 795 5394 <tel:%2B44%20%280%29151%20795%205394>,
>     Internal ext: 55394
>     M: +44 (0) 7511387895 <tel:%2B44%20%280%29%207511387895>
>     Email id: smitra at liverpool.ac.uk <mailto:smitra at liverpool.ac.uk>
>     Alternative Email id: suparna.mitra.sm at gmail.com
>     <mailto:suparna.mitra.sm at gmail.com>
>
>
>
>
> -- 
> Dr. Suparna Mitra
> Wolfson Centre for Personalised Medicine
> Department of Molecular and Clinical Pharmacology
> Institute of Translational Medicine University of Liverpool
> Block A: Waterhouse Buildings,  L69 3GL Liverpool
>
> Tel.  +44 (0)151 795 5394, Internal ext: 55394
> M: +44 (0) 7511387895
> Email id: smitra at liverpool.ac.uk <mailto:smitra at liverpool.ac.uk>
> Alternative Email id: suparna.mitra.sm at gmail.com 
> <mailto:suparna.mitra.sm at gmail.com>
>

-- 
James W. MacDonald, M.S.
Biostatistician
University of Washington
Environmental and Occupational Health Sciences
4225 Roosevelt Way NE, # 100
Seattle WA 98105-6099



More information about the Bioconductor mailing list