Hello,

I have two questions with regard to the AgiMicroRna package and they 
are FilterMicroRna and SignificantMicroRna functions. I was able to 
figure most of the package easily except these two functions. If any 
one can give me some hints or suggestions to get these two functions 
working, I will really appreciate your help. I think that the 
alternative way I used to filter control probes is causing trouble 
with the SignficantMicroRna function.

I used this package to analyze Agilent Zebrafish miRNA arrays which 
were scanned using GenePix 4000B scanner and data extracted using 
Feature Extraction Software 9.5.3. I have posted the session info 
below (two functions that need attention are hightlighted).

Session Info

 >  library("AgiMicroRna")
 >  targets.micro=readTargets(infile="targets.txt", verbose=TRUE)

Target File
                FileName Treatment GErep Subject
36_DMSO_1 36_DMSO_1.txt    36DMSO     1       1
36_DMSO_2 36_DMSO_2.txt    36DMSO     1       2
36_DMSO_3 36_DMSO_3.txt    36DMSO     1       3
36_TCDD_1 36_TCDD_1.txt    36TCDD     2       1
36_TCDD_2 36_TCDD_2.txt    36TCDD     2       2
36_TCDD_3 36_TCDD_3.txt    36TCDD     2       3
60_DMSO_1 60_DMSO_1.txt    60DMSO     3       1
60_DMSO_2 60_DMSO_2.txt    60DMSO     3       2
60_DMSO_3 60_DMSO_3.txt    60DMSO     3       3
60_TCDD_1 60_TCDD_1.txt    60TCDD     4       1
60_TCDD_2 60_TCDD_2.txt    60TCDD     4       2
60_TCDD_3 60_TCDD_3.txt    60TCDD     4       3

 > dd.micro=read.maimages(targets.micro$FileName, 
columns=list(R="gTotalGeneSignal",G= 
"gTotalProbeSignal",Rb="gMeanSignal", Gb="gProcessedSignal"), 
annotation=c("ProbeUID","ControlType","ProbeName","GeneName","SystematicName", 
"sequence", "accessions","probe_mappings", 
"gIsGeneDetected","gIsSaturated","gIsFeatNonUnifOL", 
"gIsFeatPopnOL","chr_coord","gBGMedianSignal","gBGUsed"))
Read 36_DMSO_1.txt
Read 36_DMSO_2.txt
Read 36_DMSO_3.txt
Read 36_TCDD_1.txt
Read 36_TCDD_2.txt
Read 36_TCDD_3.txt
Read 60_DMSO_1.txt
Read 60_DMSO_2.txt
Read 60_DMSO_3.txt
Read 60_TCDD_1.txt
Read 60_TCDD_2.txt
Read 60_TCDD_3.txt
 > cvArray(dd.micro, "MeanSignal", targets.micro, verbose=TRUE)
Foreground: MeanSignal

         FILTERING BY ControlType FLAG

  RAW DATA:                       5335
  PROBES without CONTROLS:        4620
----------------------------------
   (Non-CTRL) Unique Probe:  490
   (Non-CTRL) Unique Genes:  231
----------------------------------
DISTRIBUTION OF REPLICATED NonControl Probes
reps
   5   6   7  10
  20  18  36 416
------------------------------------------------------
Replication at Probe level- MEDIAN  CV
36_DMSO_1 36_DMSO_2 36_DMSO_3 36_TCDD_1 36_TCDD_2 36_TCDD_3 60_DMSO_1 
60_DMSO_2 60_DMSO_3
     0.078     0.081     0.091     0.081     0.077     0.067 
0.076     0.066     0.103
60_TCDD_1 60_TCDD_2 60_TCDD_3
     0.073     0.086     0.069
------------------------------------------------------
DISTRIBUTION OF REPLICATED Noncontrol Genes
reps
  20
231
------------------------------------------------------
 > ddTGS.rma = rmaMicroRna(dd.micro, normalize=TRUE, background=FALSE)
Calculating Expression
 > ddPROC = filterMicroRna(ddTGS.rma, dd.micro, control = TRUE, 
IsGeneDetected = TRUE, wellaboveNEG = FALSE, limIsGeneDetected = 50, 
limNEG = 25, makePLOT = FALSE, targets.micro, verbose = TRUE)
FILTERING PROBES BY FLAGS


FILTERING BY ControlType
Error in matrix(ddFILT$other$gIsGeneDetected, nrow = dim(ddFILT)[1], 
ncol = dim(ddFILT)[2]) :
   attempt to set an attribute on NULL

 > MMM = ddTGS.rma$Rb
 > colnames(MMM) = colnames(dd.micro$Rb)
 > maintitle='TGS.rma'
 > colorfill='blue'
 > ddaux=ddTGS.rma
 > ddaux$G=MMM
 > mvaMicroRna(ddaux, maintitle, verbose=TRUE)

------------------------------------------------------
mvaMicroRna info:
  FEATURES :      231
  POSITIVE CTRL:          12
  NEGATIVE CTRL:          7
  STRUCTURAL:             3
 > rm(ddaux)
 > RleMicroRna(MMM,"RLE TGS.rma", colorfill)
 > boxplotMicroRna(MMM, maintitle, colorfill)
 > plotDensityMicroRna(MMM, maintitle)
 > spottypes = readSpotTypes()
 > ddTGS.rma$genes$Status = controlStatus(spottypes, ddTGS.rma)
Matching patterns for: ProbeName GeneName
Found 231 gene
Found 1 BLANK
Found 1 Blank
Found 0 blank
Found 6 positive
Found 0 negative
Found 0 flag1
Found 0 flag2
Found 6 flag3
Found 5 flag4
Found 1 flag5
Setting attributes: values
 > i = ddTGS.rma$genes$Status == "gene"
 > esetPROC = esetMicroRna(ddTGS.rma[i,], targets.micro, 
makePLOT=TRUE, verbose = TRUE)
outPUT DATA: esetPROC
Features  Samples
      231       12
 > design=model.matrix(~-1+treatment)
 > print(design)
    treatment36DMSO treatment36TCDD treatment60DMSO treatment60TCDD
1                1               0               0               0
2                1               0               0               0
3                1               0               0               0
4                0               1               0               0
5                0               1               0               0
6                0               1               0               0
7                0               0               1               0
8                0               0               1               0
9                0               0               1               0
10               0               0               0               1
11               0               0               0               1
12               0               0               0               1
attr(,"assign")
[1] 1 1 1 1
attr(,"contrasts")
attr(,"contrasts")$treatment
[1] "contr.treatment"

 > fit=lmFit(esetPROC, design)
 > cont.matrix = makeContrasts(treatment36TCDDvstreatment36DMSO = 
treatment36TCDD-treatment36DMSO, treatment60TCDDvstreatment60DMSO = 
treatment60TCDD-treatment60DMSO,treatment60TCDDvstreatment36TCDD = 
treatment60TCDD-treatment36TCDD, treatment60DMSOvstreatment36DMSO = 
treatment60DMSO-treatment36DMSO, levels=design)
 > print(cont.matrix)
                  Contrasts
Levels            treatment36TCDDvstreatment36DMSO 
treatment60TCDDvstreatment60DMSO
   treatment36DMSO                               -1 
              0
   treatment36TCDD                                1 
              0
   treatment60DMSO                                0 
             -1
   treatment60TCDD                                0 
              1
                  Contrasts
Levels            treatment60TCDDvstreatment36TCDD 
treatment60DMSOvstreatment36DMSO
   treatment36DMSO                                0 
             -1
   treatment36TCDD                               -1 
              0
   treatment60DMSO                                0 
              1
   treatment60TCDD                                1 
              0
 > fit2 = contrasts.fit(fit,cont.matrix)
 > print(head(fit2$coeff))
             Contrasts
              treatment36TCDDvstreatment36DMSO treatment60TCDDvstreatment60DMSO
   dre-let-7a                      0.038640984                      0.013333873
   dre-let-7b                      0.074038749                     -0.031608286
   dre-let-7c                      0.026244357                     -0.005682488
   dre-let-7d                      0.067340768                      0.055567054
   dre-let-7e                      0.004569306                      0.136348664
   dre-let-7f                      0.042880109                      0.085568058
             Contrasts
              treatment60TCDDvstreatment36TCDD treatment60DMSOvstreatment36DMSO
   dre-let-7a                        1.7358343                       1.76114142
   dre-let-7b                        0.1366920                       0.24233899
   dre-let-7c                        0.9920976                       1.02402449
   dre-let-7d                        0.8098432                       0.82161694
   dre-let-7e                        0.1186829                      -0.01309647
   dre-let-7f                        1.1245878                       1.08189990
 > fit2 = eBayes(fit2)
 > fit2 = basicLimma(esetPROC, design, cont.matrix, verbose = TRUE)
DATA
Features  Samples
      231       12

 > DE = getDecideTests(fit2, DEmethod = "separate", MTestmethod = 
"BH", PVcut = 0.1, verbose = TRUE)

------------------------------------------------------
  Method for Selecting DEGs: separate
  Multiple Testing  method:  BH - pval 0.1

      treatment36TCDDvstreatment36DMSO treatment60TCDDvstreatment60DMSO
UP                                  0                                5
DOWN                                0                                1
      treatment60TCDDvstreatment36TCDD treatment60DMSOvstreatment36DMSO
UP                                 56                               51
DOWN                               80                               91
------------------------------------------------------
 > pvalHistogram(fit2, DE, PVcut = 0.1, DEmethod ="separate", 
MTestmethod="BH",cont.matrix, verbose= TRUE)
 > significantMicroRna(esetPROC, ddTGS.rma, targets.micro, fit2, 
cont.matrix, DE, DEmethod = "separate", MTestmethod= "BH", PVcut = 
0.1, Mcut=0, verbose=TRUE)
------------------------------------------------------
CONTRAST:  1  -  treatment36TCDDvstreatment36DMSO

Error in data.frame(PROBE_ID, as.character(GENE_ID), 
as.character(chr_coord),  :
   arguments imply differing number of rows: 231, 0

Thank you very much,

Sincerely, Neel


Neel Aluru Ph.D.
Post doctoral Scholar
Biology Department
Redfield 304 (MS#32)
Woods Hole Oceanographic Institution
Woods Hole MA 02543 USA
Phone: (508) 289-3607 [Office]
774-392-3727 [Cell]
RID: A-7237-2009  
	[[alternative HTML version deleted]]

