[BioC] Fwd: maSigPro - NA values for coefficient estimates - Is a polynomical fit recommended for my data?
Ana Conesa
aconesa at cipf.es
Tue Nov 10 16:20:55 CET 2009
Dear Jeremy
My apologies for the late response.
First of all: maSigPro compares profiles. This means that you must have
overlapping time points, not necessarily the same time points, but
overlapping at least. This means that if you have a group with only one time
point, you cannot do much with this. For the other 2 groups you can compare
only the segment corresponding to time points 0 and 5. You can add the rest
of the time points of the last group, but this will not be compared.
In you case I would do the following:
- Compare groups 2 and 3 considering time points 0 and 5 and using degree = 1
- Compare all groups only at time point 0 using limma
- Study the profiles changes ONLY for group 3 using a degree = 2 or 3 (if you
do not have replicates, surely not more than 2)
This, of course if this analysis makes sense for your data, but I guess that
if you designed it like this, your questions could be answered with these
analyzes.
Could it be?
Hope this answer your question. If not, come back to me..
Best regards
Ana
On Monday 09 November 2009 19:30:58 you wrote:
> Hi dear Ana,
>
> Please see my question below. I posted this last week on the BioC list
> but did not get any answer. I would appreciate if you show me a path
> here.
> I am sorry for emailing you directly. Thought I should as no body is
> replying to this question.
>
>
> ---------- Forwarded message ----------
> From: jeremy wilson <jeremy.wilson88 at gmail.com>
> Date: Thu, Nov 5, 2009 at 11:05 AM
> Subject: maSigPro - NA values for coefficient estimates - Is a
> polynomical fit recommended for my data?
> To: bioconductor at stat.math.ethz.ch
>
>
> Dear BioConductors,
>
> I have unbalanced number of time points in my experiment. In 3 groups
> I have, one group has only observations at one time point (time 0),
> the other at two time points (times 0, 5) and the last group with four
> time points(0, 5, 12, 31). I am wondering if I can use maSigPro for
> this type of data. In this case, what order of polynomial fit is good?
> Do you recommend to use degree = 1 to reduce the degree of the model
> as the best model is the one with the least order? or is it good to
> model with degree = 3? I have 4 time points in one group.
>
> Is it better to do this kind of unbalanced time points data with
> LIMMA? There are many comparisons to consider if I use LIMMA. What do
> you suggest to use?
>
> I did fit a degree = 3 polynomial model and got NA values. The NA
> values appeared as below. There were no NA values for the remaining
> coefficients. Are these due to no values for some experimental groups
> at some time points? Can I ignore these warnings?
>
> Here are my commands and warnings I got..
> design.init<-read.table("maSigPro_PhenotypeInput.txt", sep="\t",
> header=TRUE, row.names=1)
>
> > design<-make.design.matrix(design.init, deg=3) # Create a maSigPro design
> > matrix
> >
> > design
>
> $dis
> VehiclevsNaive BCGvsNaive Time TimexVehicle TimexBCG
> Time2 Time2xVehicle Time2xBCG Time3 Time3xVehicle Time3xBCG
> Veh Day0 A2.CEL 1 0 0 0 0
> 0 0 0 0 0 0
> Veh Day0 A4.CEL 1 0 0 0 0
> 0 0 0 0 0 0
> Veh Day0 A5.CEL 1 0 0 0 0
> 0 0 0 0 0 0
> Veh Day12 A1.CEL 1 0 12 12 0
> 144 144 0 1728 1728 0
> Veh Day12 A3.CEL 1 0 12 12 0
> 144 144 0 1728 1728 0
> Veh Day0 B1.CEL 1 0 0 0 0
> 0 0 0 0 0 0
> Veh Day0 B3.CEL 1 0 0 0 0
> 0 0 0 0 0 0
> Veh Day0 B4.CEL 1 0 0 0 0
> 0 0 0 0 0 0
> Veh Day12 B2.CEL 1 0 12 12 0
> 144 144 0 1728 1728 0
> Veh Day12 B3.CEL 1 0 12 12 0
> 144 144 0 1728 1728 0
> Veh Day0 C1.CEL 1 0 0 0 0
> 0 0 0 0 0 0
> Veh Day0 C3.CEL 1 0 0 0 0
> 0 0 0 0 0 0
> Veh Day0 C4.CEL 1 0 0 0 0
> 0 0 0 0 0 0
> Veh Day0 C5.CEL 1 0 0 0 0
> 0 0 0 0 0 0
> Veh Day12 C1.CEL 1 0 12 12 0
> 144 144 0 1728 1728 0
> Veh Day12 C3.CEL 1 0 12 12 0
> 144 144 0 1728 1728 0
> Day 0 A1.CEL 0 1 0 0 0
> 0 0 0 0 0 0
> Day 0 A2.CEL 0 1 0 0 0
> 0 0 0 0 0 0
> Day 0 A3.CEL 0 1 0 0 0
> 0 0 0 0 0 0
> Day 5 A1.CEL 0 1 5 0 5
> 25 0 25 125 0 125
> Day 5 A2.CEL 0 1 5 0 5
> 25 0 25 125 0 125
> Day 5 A3.CEL 0 1 5 0 5
> 25 0 25 125 0 125
> Day 12 A1.CEL 0 1 12 0 12
> 144 0 144 1728 0 1728
> Day 12 A2.CEL 0 1 12 0 12
> 144 0 144 1728 0 1728
> Day 12 A3.CEL 0 1 12 0 12
> 144 0 144 1728 0 1728
> Day 31 A1.CEL 0 1 31 0 31
> 961 0 961 29791 0 29791
> Day 31 A2.CEL 0 1 31 0 31
> 961 0 961 29791 0 29791
> Day 31 A3.CEL 0 1 31 0 31
> 961 0 961 29791 0 29791
> Day 0 B1.CEL 0 1 0 0 0
> 0 0 0 0 0 0
> Day 0 B2.CEL 0 1 0 0 0
> 0 0 0 0 0 0
> Day 0 B3.CEL 0 1 0 0 0
> 0 0 0 0 0 0
> Day 5 B1.CEL 0 1 5 0 5
> 25 0 25 125 0 125
> Day 5 B2.CEL 0 1 5 0 5
> 25 0 25 125 0 125
> Day 5 B3.CEL 0 1 5 0 5
> 25 0 25 125 0 125
> Day 12 B1.CEL 0 1 12 0 12
> 144 0 144 1728 0 1728
> Day 12 B2.CEL 0 1 12 0 12
> 144 0 144 1728 0 1728
> Day 12 B3.CEL 0 1 12 0 12
> 144 0 144 1728 0 1728
> Day 31 B1.CEL 0 1 31 0 31
> 961 0 961 29791 0 29791
> Day 31 B2.CEL 0 1 31 0 31
> 961 0 961 29791 0 29791
> Day 31 B3.CEL 0 1 31 0 31
> 961 0 961 29791 0 29791
> Day 0 C1.CEL 0 1 0 0 0
> 0 0 0 0 0 0
> Day 0 C2.CEL 0 1 0 0 0
> 0 0 0 0 0 0
> Day 0 C3.CEL 0 1 0 0 0
> 0 0 0 0 0 0
> Day 0 C4.CEL 0 1 0 0 0
> 0 0 0 0 0 0
> Day 5 C1.CEL 0 1 5 0 5
> 25 0 25 125 0 125
> Day 5 C2.CEL 0 1 5 0 5
> 25 0 25 125 0 125
> Day 5 C3.CEL 0 1 5 0 5
> 25 0 25 125 0 125
> Day 12 C1.CEL 0 1 12 0 12
> 144 0 144 1728 0 1728
> Day 12 C2.CEL 0 1 12 0 12
> 144 0 144 1728 0 1728
> Day 12 C3.CEL 0 1 12 0 12
> 144 0 144 1728 0 1728
> Day 31 C1.CEL 0 1 31 0 31
> 961 0 961 29791 0 29791
> Day 31 C2.CEL 0 1 31 0 31
> 961 0 961 29791 0 29791
> Day 31 C3.CEL 0 1 31 0 31
> 961 0 961 29791 0 29791
> Naive A1.CEL 0 0 0 0 0
> 0 0 0 0 0 0
> Naive A2.CEL 0 0 0 0 0
> 0 0 0 0 0 0
> Naive B1.CEL 0 0 0 0 0
> 0 0 0 0 0 0
> Naive B2.CEL 0 0 0 0 0
> 0 0 0 0 0 0
> Naive C1.CEL 0 0 0 0 0
> 0 0 0 0 0 0
> Naive C2.CEL 0 0 0 0 0
> 0 0 0 0 0 0
> Naive C3.CEL 0 0 0 0 0
> 0 0 0 0 0 0
>
> $groups.vector
> [1] "VehiclevsNaive" "BCGvsNaive" "Naive"
> "VehiclevsNaive" "BCGvsNaive" "Naive" "VehiclevsNaive"
> "BCGvsNaive"
> [9] "Naive" "VehiclevsNaive" "BCGvsNaive"
>
> $edesign
> Time Replicate Naive Vehicle BCG
> Veh Day0 A2.CEL 0 1 0 1 0
> Veh Day0 A4.CEL 0 1 0 1 0
> Veh Day0 A5.CEL 0 1 0 1 0
> Veh Day12 A1.CEL 12 2 0 1 0
> Veh Day12 A3.CEL 12 2 0 1 0
> Veh Day0 B1.CEL 0 3 0 1 0
> Veh Day0 B3.CEL 0 3 0 1 0
> Veh Day0 B4.CEL 0 3 0 1 0
> Veh Day12 B2.CEL 12 4 0 1 0
> Veh Day12 B3.CEL 12 4 0 1 0
> Veh Day0 C1.CEL 0 5 0 1 0
> Veh Day0 C3.CEL 0 5 0 1 0
> Veh Day0 C4.CEL 0 5 0 1 0
> Veh Day0 C5.CEL 0 5 0 1 0
> Veh Day12 C1.CEL 12 6 0 1 0
> Veh Day12 C3.CEL 12 6 0 1 0
> Day 0 A1.CEL 0 7 0 0 1
> Day 0 A2.CEL 0 7 0 0 1
> Day 0 A3.CEL 0 7 0 0 1
> Day 5 A1.CEL 5 8 0 0 1
> Day 5 A2.CEL 5 8 0 0 1
> Day 5 A3.CEL 5 8 0 0 1
> Day 12 A1.CEL 12 9 0 0 1
> Day 12 A2.CEL 12 9 0 0 1
> Day 12 A3.CEL 12 9 0 0 1
> Day 31 A1.CEL 31 10 0 0 1
> Day 31 A2.CEL 31 10 0 0 1
> Day 31 A3.CEL 31 10 0 0 1
> Day 0 B1.CEL 0 11 0 0 1
> Day 0 B2.CEL 0 11 0 0 1
> Day 0 B3.CEL 0 11 0 0 1
> Day 5 B1.CEL 5 12 0 0 1
> Day 5 B2.CEL 5 12 0 0 1
> Day 5 B3.CEL 5 12 0 0 1
> Day 12 B1.CEL 12 13 0 0 1
> Day 12 B2.CEL 12 13 0 0 1
> Day 12 B3.CEL 12 13 0 0 1
> Day 31 B1.CEL 31 14 0 0 1
> Day 31 B2.CEL 31 14 0 0 1
> Day 31 B3.CEL 31 14 0 0 1
> Day 0 C1.CEL 0 15 0 0 1
> Day 0 C2.CEL 0 15 0 0 1
> Day 0 C3.CEL 0 15 0 0 1
> Day 0 C4.CEL 0 15 0 0 1
> Day 5 C1.CEL 5 16 0 0 1
> Day 5 C2.CEL 5 16 0 0 1
> Day 5 C3.CEL 5 16 0 0 1
> Day 12 C1.CEL 12 17 0 0 1
> Day 12 C2.CEL 12 17 0 0 1
> Day 12 C3.CEL 12 17 0 0 1
> Day 31 C1.CEL 31 18 0 0 1
> Day 31 C2.CEL 31 18 0 0 1
> Day 31 C3.CEL 31 18 0 0 1
> Naive A1.CEL 0 19 1 0 0
> Naive A2.CEL 0 19 1 0 0
> Naive B1.CEL 0 20 1 0 0
> Naive B2.CEL 0 20 1 0 0
> Naive C1.CEL 0 21 1 0 0
> Naive C2.CEL 0 21 1 0 0
> Naive C3.CEL 0 21 1 0 0
>
> > fit <- p.vector(efiltered.mat, design, Q = 0.05, MT.adjust = "BH",min.obs
> > = 20) tstep <- T.fit(fit, step.method = "two.ways.backward", alfa = 0.05)
>
> I got the following warnings after the above step
> warnings:
> 1: In rbind(sol, as.numeric(c(p.value, result$r.squared, ... :
> NAs introduced by coercion
> 2: In rbind(sol, as.numeric(c(p.value, result$r.squared, ... :
> NAs introduced by coercion
>
> >sigs <- get.siggenes(tstep, rsq = 0.6, vars = "groups")
> > see.genes(sigs$sig.genes$BCGvsNaive, main = "BCGvsNaive", show.fit =
> > T,dis =design$dis, cluster.method="kmeans" ,cluster.data = 1, k = 9)
>
> Warning messages:
> 1: In rbind(sol, as.numeric(c(p.value, result$r.squared, p.valor))) :
> NAs introduced by coercion
> 2: In rbind(sol, as.numeric(c(p.value, result$r.squared, p.valor))) :
> NAs introduced by coercion
>
>
>
> betaTimexBCG betaTime2 betaTime2xVehicle
> betaTime2xBCG betaTime3
> AB002558_at NA 0.0000000000 NA
> NA -1.824498e-05
> AF001953_g_at NA 0.0114056076 NA
> NA -2.616775e-04
> AF067795_at NA 0.0011159473 NA
> NA -3.539580e-05
> D13120_s_at NA 0.0011509775 NA
> NA 0.000000e+00
> E01534cds_s_at NA -0.0029700508 NA
> NA 8.818569e-05
> L24896_s_at NA -0.0164766816 NA
> NA 3.728959e-04
> M20156_at NA -0.0011509235 NA
> NA 0.000000e+00
> M25638_s_at NA 0.0137020816 NA
> NA -3.174337e-04
> M58404_at NA -0.0123769559 NA
> NA 2.928362e-04
> M60921_at NA 0.0000000000 NA
> NA -1.465997e-05
> M91652complete_seq_at NA -0.0020565596 NA
> NA 0.000000e+00
> M94918mRNA_f_at NA -0.0256353952 NA
> NA 6.500460e-04
> M94919mRNA_f_at NA 0.0000000000 NA
> NA 3.851441e-05
> rc_AA859783_at NA 0.0009030979 NA
> NA 0.000000e+00
> rc_AA875171_g_at NA -0.0007324529 NA
> NA 0.000000e+00
> rc_AA891171_s_at NA 0.0010095092 NA
> NA 0.000000e+00
> rc_AA891727_g_at NA -0.0047867258 NA
> NA 1.441130e-04
> rc_AA892123_at NA -0.0038206649 NA
> NA 1.165693e-04
> rc_AA892248_g_at NA 0.0000000000 NA
> NA 0.000000e+00
> rc_AA893246_at NA 0.0010566522 NA
> NA 0.000000e+00
> rc_AA924772_at NA -0.0008263709 NA
> NA 0.000000e+00
> rc_AA944324_at NA -0.0013202929 NA
> NA 0.000000e+00
> rc_AI007824_g_at NA 0.0000000000 NA
> NA 0.000000e+00
> rc_AI013194_at NA 0.0000000000 NA
> NA 0.000000e+00
> rc_AI102044_at NA 0.0000000000 NA
> NA 0.000000e+00
> rc_AI104544_at NA 0.0019113263 NA
> NA 0.000000e+00
> rc_AI170268_at NA -0.0012527673 NA
> NA 0.000000e+00
> rc_AI176662_s_at NA 0.0036907561 NA
> NA -1.149829e-04
> rc_AI179576_s_at NA -0.0275651514 NA
> NA 6.926836e-04
> rc_AI235364_at NA -0.0028867693 NA
> NA 8.997026e-05
> rc_AI237836_at NA 0.0000000000 NA
> NA 0.000000e+00
> S68135_s_at NA 0.0064503716 NA
> NA -1.943265e-04
> S74351_s_at NA 0.0066852568 NA
> NA -2.153613e-04
> U02553cds_s_at NA 0.0000000000 NA
> NA -2.903164e-05
> U12568_at NA -0.0012220408 NA
> NA 0.000000e+00
> U19866_at NA -0.0010647292 NA
> NA 0.000000e+00
> U39875_at NA -0.0007813637 NA
> NA 0.000000e+00
> U84727_at NA 0.0104752791 NA
> NA -2.460059e-04
> X03347cds_g_at NA 0.0022984842 NA
> NA -7.116160e-05
> X54510_g_at NA 0.0000000000 NA
> NA 2.572912e-05
> X56325mRNA_s_at NA -0.0214234307 NA
> NA 5.057676e-04
> X58389cds_s_at NA -0.0034907821 NA
> NA 1.055627e-04
> X63375exon_at NA 0.0132059943 NA
> NA -3.082255e-04
> X66369_at NA 0.0014723516 NA
> NA 0.000000e+00
> X94242_at NA 0.0000000000 NA
> NA 0.000000e+00
> betaTime3xVehicle betaTime3xBCG
> AB002558_at NA NA
> AF001953_g_at NA NA
> AF067795_at NA NA
> D13120_s_at NA NA
> E01534cds_s_at NA NA
> L24896_s_at NA NA
> M20156_at NA NA
> M25638_s_at NA NA
> M58404_at NA NA
> M60921_at NA NA
> M91652complete_seq_at NA NA
> M94918mRNA_f_at NA NA
> M94919mRNA_f_at NA NA
> rc_AA859783_at NA NA
> rc_AA875171_g_at NA NA
> rc_AA891171_s_at NA NA
> rc_AA891727_g_at NA NA
> rc_AA892123_at NA NA
> rc_AA892248_g_at NA NA
> rc_AA893246_at NA NA
> rc_AA924772_at NA NA
> rc_AA944324_at NA NA
> rc_AI007824_g_at NA NA
> rc_AI013194_at NA NA
> rc_AI102044_at NA NA
> rc_AI104544_at NA NA
> rc_AI170268_at NA NA
> rc_AI176662_s_at NA NA
> rc_AI179576_s_at NA NA
> rc_AI235364_at NA NA
> rc_AI237836_at NA NA
> S68135_s_at NA NA
> S74351_s_at NA NA
> U02553cds_s_at NA NA
> U12568_at NA NA
> U19866_at NA NA
> U39875_at NA NA
> U84727_at NA NA
> X03347cds_g_at NA NA
> X54510_g_at NA NA
> X56325mRNA_s_at NA NA
> X58389cds_s_at NA NA
> X63375exon_at NA NA
> X66369_at NA NA
> X94242_at NA NA
>
>
> The following coefficients had no "NA" values.
> tstep$coefficients$beta0
> tstep$coefficients$betaVehiclevsNaive
> tstep$coefficients$betaBCGvsNaive
> tstep$coefficients$betaTime tstep$coefficients$betaTimexVehicle
>
> > sessionInfo()
>
> R version 2.10.0 (2009-10-26)
> i386-pc-mingw32
>
> locale:
> [1] LC_COLLATE=English_United States.1252 LC_CTYPE=English_United
> States.1252
> [3] LC_MONETARY=English_United States.1252 LC_NUMERIC=C
> [5] LC_TIME=English_United States.1252
>
> attached base packages:
> [1] tcltk stats graphics grDevices utils datasets
> methods base
>
> other attached packages:
> [1] maSigPro_1.17.5 DynDoc_1.23.1
> widgetTools_1.23.1 genefilter_1.26.4
> [5] arrayQualityMetrics_2.4.1 affyPLM_1.22.0
> preprocessCore_1.7.9 gcrma_2.17.4
> [9] rgu34acdf_2.5.0 affy_1.23.12 Biobase_2.5.8
>
> loaded via a namespace (and not attached):
> [1] affyio_1.13.5 annotate_1.23.4 AnnotationDbi_1.7.20
> beadarray_1.13.9
> [5] Biostrings_2.13.54 DBI_0.2-4 grid_2.10.0
> hwriter_1.1
> [9] IRanges_1.3.99 KernSmooth_2.23-3 lattice_0.17-26
> latticeExtra_0.6-3
> [13] limma_3.0.3 marray_1.23.0 Mfuzz_2.3.4
> RColorBrewer_1.0-2
> [17] RSQLite_0.7-3 simpleaffy_2.21.3 splines_2.10.0
> stats4_2.10.0
> [21] survival_2.35-7 tkWidgets_1.23.2 tools_2.10.0
> vsn_3.13.7
> [25] xtable_1.5-5
>
>
> I would really appreciate for any suggestions. Thanks in advance..
--
Ana Conesa
Bioinformatics and Genomics Department
Centro de Investigaciones Principe Felipe
Avda. Autopista Saler 16,
46012 Valencia, Spain
Phone: +34 96 328 96 80
Fax: +34 96 328 97 01
http://bioinfo.cipf.es/aconesa
http://www.blast2go.org
==========================================
FIRST INTERNATIONAL COURSE IN AUTOMATED
FUNCTIONAL ANNOTATION AND DATA MINING
Valencia/Orlando, September/October 2009
http://bioinfo.cipf.es/blast2gocourse
==========================================
More information about the Bioconductor
mailing list