On Fri, May 14, 2010 at 1:47 PM, Thomas Levine <thomas.lev...@gmail.com> wrote: > Creating the 5 indicator variables will be easy if you post your code > and sample data. This may also allow people to help with the first > problem you were having.
Here you go. Fragment of data file with 4000 or so pts: RawTime Time OxpcA THCA OxyA DeOxyA Aux1 Aux2 Aux3 9213290 0 71.2 51.47 36.64 14.82 3.88 3.88 0.35 9213390 0.1 66.1 46.26 30.59 15.67 3.88 3.88 0.35 9213490 0.2 81.8 54.07 44.25 9.815 3.88 3.88 0.35 9213590 0.3 58.1 56.7 32.96 23.74 3.88 3.88 0.35 9213690 0.4 50.2 48.55 24.39 24.16 3.88 3.88 0.35 9213790 0.5 59.7 47.02 28.08 18.94 3.88 3.88 0.35 9213890 0.6 68.0 52.78 35.9 16.88 3.88 3.88 0.35 9213990 0.7 60.7 46.62 28.28 18.34 3.88 3.88 0.35 9214090 0.8 57.5 47.81 27.51 20.31 3.88 3.88 0.35 9214190 0.9 60.4 52.13 31.49 20.64 3.88 3.88 0.35 9214290 1 25.5 48.21 12.31 35.9 3.88 3.88 0.35 9214390 1.1 58.3 49.67 28.95 20.72 3.88 3.88 0.35 9214490 1.2 65.5 57.28 37.52 19.76 3.88 3.88 0.35 9214590 1.3 75.1 53.84 40.44 13.4 3.88 3.88 0.35 9214690 1.4 58.9 40.41 23.79 16.62 3.88 3.88 0.35 9214790 1.5 76.6 50.1 38.38 11.73 3.88 3.88 0.35 9214890 1.6 72.4 52.78 38.21 14.57 3.88 3.88 0.35 9214990 1.7 66.5 41.27 27.43 13.84 3.88 3.88 0.35 9215090 1.8 62.4 52.4 32.68 19.72 3.88 3.88 0.35 9215190 1.9 43.9 45.11 19.8 25.31 3.88 3.88 0.35 9215290 2 47.1 46.17 21.75 24.42 3.88 3.88 0.35 #convert binary Aux1--Aux3 into decimal x<- (d$Aux1>mean(range(d$Aux1)))*1 + (d$Aux2>mean(range(d$Aux2)))*2 +(d$Aux3>mean(range(d$Aux3)))*4 x<-abs(x-3) #do this because resp to +contrast is same as to -contrast hbr<-d$DeOxyA #detrend() is my own function which detrends :-) hbr<-detrend(hbr) #lm.lag is my own function which shifts x in time so it produces the best fit with y #then fits lm(y~(shifted x)) #This treats x as continuous (which it is) fit1<-lm.lag(x,hbr,lag=c(50,200)); summary(fit1$fit); fit1$lag xs<-shift(x,lag=97) plot(0:3, aggregate(cbind(hbr,xs), list(xs), mean)[,2]); abline(a=18.00991,b= -.18183) #this is the simple indicator var approach that compares 0 vs 1, 0 vs 2, 0 vs 3 # It is NOT what I want fit2<-lm(hbr~factor(xs),contrasts=contr.treatment(3)) summary(fit2) I am wondering, first of all, if an lm() approach can tell me what I want to know. Namely, the coefficients of the linear fit for y~x, and also some coefficients which indicate what kind of nonlinearities I have and how strong they are. I was asking how to do it and threw out lm() as a possibility. I am open to anything. Thanks for any help Bill [[alternative HTML version deleted]] ______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.