[BioC] can't get summary.stat after "invariantset" normalizat ion

[Laurent Gautier]

Well... seems I've been boasting too much (and I was secretely hoping
it was the devel version...). I do not have run affy-release on 1.6.0,
I can only suggests reasons.

> On Thu, Oct 17, 2002 at 11:24:24AM -0400, Man, Michael wrote:
> > Any idea?  BTW, it works fine when I used unnormalized data.  -michael

This is the funny part. I wonder why it does 
> 
> 
> Several even... which version of R and of the package are you using ?
> 
> 
> L.
> 
> 
> > 
> > ###############################
> > > tmp <- normalize.Plob.invariantset(latin)  # use the first CHIP as
> > reference

The following do it too (preferred way):

tmp <- normalize(latin, method="invariantset")


> > > latin.LiWong <- express(tmp, normalize=F, bg=subtractmm,
> > summary.stat=li.wong)
> > Background correcting
> > Preparing Data
> > Computing expression. This may take a while.
> > Error in while (change > delta & iter < maxit) { : 
> >         missing value where logical needed

Would you have variables called 'T' or 'F' around ?
(the code relied (too much) on (T == TRUE) and (F == FALSE) but this
was fixed about a week ago I believe). 

> > In addition: Warning message: 
> > No convergence in inner loop after 50 in outerler tieration 4 
> >  in: fit.li.wong(t(data.matrix), remove.outliers, normal.array.quantile,  
> > 

This only a warning.

> > 
> > > latin.LiWong <- express(tmp, normalize=F, summary.stat=li.wong)
> > Background correcting
> > Error in density(x, kernel = "epanechnikov", n = n.pts) : 
> >         x contains missing values
> >

In this case the default 'bg' is rma. You might get 'masked' probes from
your CEL files (to avoid taking them into consideration, try to set the
flags related to masked and outlying probes to FALSE (see man related
man page)). I thought it was fixed... but may be only in affy-devel.
Try to patch bg.parameters with the following:

bg.parameters <- function(pm, mm, n.pts=2^14){

  max.density <- function(x,n.pts){
    aux <- density(x, kernel="epanechnikov", n=n.pts, na.rm=TRUE) 
    aux$x[order(-aux$y)[1]] 
  }
    mmbg <- max.density(mm,n.pts)
    pmbg <- max.density(pm,n.pts)
		   
    bg.data <- mm[mm < mmbg]
    bg.data <- bg.data-mmbg
    bgsd <- sqrt(sum(bg.data^2)/(length(bg.data)-1))*sqrt(2)/.85
			  
    sig.data <- pm[pm > pmbg]
    sig.data <- sig.data-pmbg
			       
    alpha <- 1/mean(sig.data)
    mubg <- mmbg
    list(alpha=alpha,mu=mubg,sigma=bgsd)  
}
				     

> > 
> > > latin.avdf <- express(tmp,normalize=F,bg=subtractmm,summary.stat=avdiff)
> > Background correcting
> > Preparing Data
> > Computing expression. This may take a while.
> > Error in var(x, na.rm = na.rm) : missing observations in cov/cor
> > 

A patch was contributed for that and I was convinced I applied it to
both the released and devel version (I cannot find track of it in CVS ?!
I am (much) confused now...)
The reason is like above: presence of NAs.

> > 
> > > latin.avdf <- express(tmp,normalize=F,bg=subtractmm,
> > +                       summary.stat=function(x) apply(x,2,mean,trim=.2))
> > Background correcting
> > Preparing Data
> > Computing expression. This may take a while.
> > Error in psort(x, partial) : index 13 outside bounds


What does give 'traceback()' on this one ?




L.