Re: [R] nls and factor
Thanks, it was actually p.249, at least in my MASS3. but that solved my doubt. I've have another doubt, can this factor interact with one of the parameters in the model? My problem is basically a Michaelis Menten term, where this factor determines a different Km. The rest of the parameters in the model are the same. But I don't know how to write the nls formula, or if it is possible. This is a toy example, B and A are the two factors, conc is the concentration and t is the temperature: ## Generate independent variables Bconc-runif(30,0.1,10) Aconc-runif(30,0.1,10) At-runif(30,1,30) Bt-runif(30,1,30) ## These are the parameters I want to calculate from ## my real data BKm-1 AKm-0 EBoth--0.41 # These are my simulated dependent variables yB-100*exp(EBoth*Bt)*Bconc/(BKm+Bconc)+rnorm(30,0,1) yA-75*exp(EBoth*At)*Aconc/(AKm+Aconc)+rnorm(30,0,1) #The separate models BModel-nls(Response~lev*exp(Ev*t)*conc/(Km+conc),data=list(Response=yB,t=Bt,conc=Bconc),start=list(lev=90,Ev=-0.5,Km=0.8),trace=TRUE) AModel-nls(Response~lev*exp(Ev*t)*conc/(Km+conc),data=list(Response=yA,t=At,conc=Aconc),start=list(lev=90,Ev=-0.5,Km=0.8),trace=TRUE) ## I want to obtain a combined model of the form: ## Y=Intercept[1:2]*exp(Eboth*t)*conc/(Km[1:2]+conc) ## where I have a common E but two intercepts and two ## Kms (one of them should in fact be zero) yBoth-c(yB,yA) concBoth-c(Bconc,Aconc) tBoth-c(At,Bt) AorB-as.factor(c(rep(0,length(yA)),rep(1,length(yB ## Amongst other things I've tried FullModel-nls(Response~lev[AorB]*exp(Ev*t)*conc/(Km[AorB]+conc),data=list(Response=yBoth,t=tBoth,conc=concBoth),start=list(lev=c(90,70),Ev=-0.5,Km=c(0.8,0)),trace=TRUE) ## but i get to a singular gradient Any other pointers, thanks Manuel --- Prof Brian Ripley [EMAIL PROTECTED] escribió: On Thu, 20 Apr 2006, Manuel Gutierrez wrote: Is it possible to include a factor in an nls formula? Yes. What do you intend by it? If you mean what it would mean for a lm formula, you need A[a] and starting values for A. There's an example on p.219 of MASS4. I've searched the help pages without any luck so I guess it is not feasible. I've given it a few attempts without luck getting the message: + not meaningful for factors in: Ops.factor(independ^EE, a) This is a toy example, my realworld case is much more complicated (and can not be solved linearizing an using lm) a-as.factor(c(rep(1,50),rep(0,50))) independ-rnorm(100) respo-rep(NA,100) respo[a==1]-(independ[a==1]^2.3)+2 respo[a==0]-(independ[a==0]^2.1)+3 nls(respo~independ^EE+a,start=list(EE=1.8),trace=TRUE) Any pointers welcomed Many Thanks, Manu __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] nls and factor
Is it possible to include a factor in an nls formula? I've searched the help pages without any luck so I guess it is not feasible. I've given it a few attempts without luck getting the message: + not meaningful for factors in: Ops.factor(independ^EE, a) This is a toy example, my realworld case is much more complicated (and can not be solved linearizing an using lm) a-as.factor(c(rep(1,50),rep(0,50))) independ-rnorm(100) respo-rep(NA,100) respo[a==1]-(independ[a==1]^2.3)+2 respo[a==0]-(independ[a==0]^2.1)+3 nls(respo~independ^EE+a,start=list(EE=1.8),trace=TRUE) Any pointers welcomed Many Thanks, Manu __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] nls and factor
Thanks Andrew. I am now trying but without much success. I don't now how to give start values for the factor?. Could you give me an example solution with my toy example? a-as.factor(c(rep(1,50),rep(0,50))) independ-1:100 respo-rep(NA,100) respo[a==1]-(independ[a==1]^2.3)+2 respo[a==0]-(independ[a==0]^2.1)+3 library(nlme) gnls(respo~independ^b+a,start=list(b=1.8)) Many thanks. Manu --- Andrew Robinson [EMAIL PROTECTED] escribió: Manuel, I don't think that it works very easily. Instead, try gnls() in the nlme package. Cheers Andrew On Thu, Apr 20, 2006 at 11:18:02AM +0200, Manuel Gutierrez wrote: Is it possible to include a factor in an nls formula? I've searched the help pages without any luck so I guess it is not feasible. I've given it a few attempts without luck getting the message: + not meaningful for factors in: Ops.factor(independ^EE, a) This is a toy example, my realworld case is much more complicated (and can not be solved linearizing an using lm) a-as.factor(c(rep(1,50),rep(0,50))) independ-rnorm(100) respo-rep(NA,100) respo[a==1]-(independ[a==1]^2.3)+2 respo[a==0]-(independ[a==0]^2.1)+3 nls(respo~independ^EE+a,start=list(EE=1.8),trace=TRUE) Any pointers welcomed Many Thanks, Manu __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Andrew Robinson Department of Mathematics and Statistics Tel: +61-3-8344-9763 University of Melbourne, VIC 3010 Australia Fax: +61-3-8344-4599 Email: [EMAIL PROTECTED] http://www.ms.unimelb.edu.au __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] se.fit in predict.nls
Sorry to be so persistent but I really need to get some measure of the error in the predictions of my nls model. I've tried to find out what predict.nls does and I've got down to MiModelo$m$predict function (newdata = list(), qr = FALSE) { eval(form[[3]], as.list(newdata), env) } environment: 0x88a076c But I can not find what is form. Any help, please. Manuel Manuel Gutierrez wrote: The option se.fit in predict.nls is currently ignored. Is there any other function available to calculate the error in the predictions? Thanks, Manuel __ LLama Gratis a cualquier PC del Mundo. Llamadas a fijos y móviles desde 1 céntimo por minuto. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] se.fit in predict.nls
The option se.fit in predict.nls is currently ignored. Is there any other function available to calculate the error in the predictions? Thanks, Manuel __ LLama Gratis a cualquier PC del Mundo. Llamadas a fijos y móviles desde 1 céntimo por minuto. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] error propagation
Are there any functions to do error propagation in R? I have done a search with little success. Any pointers to read about this topic are greatly welcomed. My specific problem is: I use a linear model (lm) to predict the biomass of an individual in a population, then I add up the biomass of all individuals to calculate the total mass of the population. I want to calculate the error in the final estimate of the total biomass. Also, I have another case where I use a nls model instead of the lm, but in nls se.fit is currently ignored as stated in the help page. Is there an alternative way to calculate the errors in the predictions from nls? Thanks, Manuel __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] addition of error terms
Are there any functions to do error propagation in R? I have done a search with little success. Any pointers to read about this topic are greatly welcomed. My specific problem is: I use a linear model (lm) to predict the biomass of an individual in a population, then I add up the biomass of all individuals to calculate the total mass of the population. I want to calculate the error in the final estimate of the total biomass. Also, I have another case where I use a nls model instead of the lm, but in nls se.fit is currently ignored as stated in the help page. Is there an alternative way to calculate the errors in the predictions from nls? Thanks, Manuel __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] dealing with multicollinearity
I have a linear model y~x1+x2 of some data where the coefficient for x1 is higher than I would have expected from theory (0.7 vs 0.88) I wondered whether this would be an artifact due to x1 and x2 being correlated despite that the variance inflation factor is not too high (1.065): I used perturbation analysis to evaluate collinearity library(perturb) P-perturb(A,pvars=c(x1,x2),prange=c(1,1)) summary(P) Perturb variables: x1 normal(0,1) x2 normal(0,1) Impact of perturbations on coefficients: mean s.d. min max (Intercept) -26.0670.270 -27.235 -25.481 x1 0.7260.0250.6720.882 x2 0.0600.0110.0370.082 I get a mean for x1 of 0.726 which is closer to what is expected. I am not an statistical expert so I'd like to know if my evaluation of the effects of collinearity is correct and in that case any solutions to obtain a reliable linear model. Thanks, Manuel Some more detailed information: A-lm(y~x1+x2) summary(A) Call: lm(formula = y ~ x1 + x2) Residuals: Min1QMedian3Q Max -4.221946 -0.484055 -0.004762 0.397508 2.542769 Coefficients: Estimate Std. Error t value Pr(|t|) (Intercept) -27.234720.27996 -97.282 2e-16 *** x10.882020.02475 35.639 2e-16 *** x20.081800.01239 6.604 2.53e-10 *** --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 Residual standard error: 0.823 on 241 degrees of freedom Multiple R-Squared: 0.8411, Adjusted R-squared: 0.8398 F-statistic: 637.8 on 2 and 241 DF, p-value: 2.2e-16 cor.test(x1,x2) Pearson's product-moment correlation data: x1 and x2 t = -3.9924, df = 242, p-value = 8.678e-05 alternative hypothesis: true correlation is not equal to 0 95 percent confidence interval: -0.3628424 -0.1269618 sample estimates: cor -0.248584 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Re: testing slopes different than a given value
Thanks to all, Yes, I meant a single test for both coefficients. Peter's reply is what I wanted. I've tried with linear.hypothesis but I must confess that with my limited statistical experience and without the car book at hand, the nomenclature for the function was a bit difficult to understand for me. A toy example of linear.hypothesis for my case would be great. Thanks, Manuel --- Peter Dalgaard [EMAIL PROTECTED] escribió: John Fox [EMAIL PROTECTED] writes: Dear Vito, Since Manuel says that he wants to obtain a test and not obtain two tests, I assume that he's interested in the F-test for the hypothesis that both coefficients are simultaneously equal to the specified values rather than in the t-tests for the individual hypotheses. Regards, John ...in which case one answer is this: y-3+0.6*x1+0.3*x2 + rnorm(100,sd=.1) # as meant, no? fm-lm(y~x1+x2) anova(fm, lm(y~offset(0.6*x1+0.3*x2))) Analysis of Variance Table Model 1: y ~ x1 + x2 Model 2: y ~ offset(0.6 * x1 + 0.3 * x2) Res.Df RSS Df Sum of Sq F Pr(F) 1 97 1.06118 2 99 1.06184 -2 -0.00066 0.0302 0.9703 -- O__ Peter Dalgaard Blegdamsvej 3 c/ /'_ --- Dept. of Biostatistics 2200 Cph. N (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907 __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] testing slopes different than a given value
In a multiple linear regression with two independent variables is there any function in R to test for the coefficients being different than some given values? Example: x1-rnorm(100) x2-rnorm(100) y-3+0.6*x1+0.3*x2 lm(y~x1+x2) Obtain a test for the coefficients for x1 being different than 0.6 and for x2 different than 0.3 Thanks Manuel __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] Rare Cases and SOM
I am trying to understand how the SOM algorithm works using library(class) SOM function. I have a 1000*10 matrix and I want to be able to summarize the different types of 10-element vectors. In my real world case it is likely that most of the 1000 values are of one kind the rest of other (this is an oversimplification). Say for example: InputA-matrix(cos(1:10),nrow=900,ncol=10,byrow=TRUE) InputB-matrix(sin(5:14),nrow=100,ncol=10,byrow=TRUE) Input-rbind(InputA,InputB) I though that a small grid of 3*3 would be enough to extract the patterns in such simple matrix : GridWidth-3 GridLength-3 gr - somgrid(xdim=GridWidth,ydim=GridLength,topo = hexagonal) test.som - SOM(Input, gr) par(mfrow=c(GridLength,GridWidth)) for(i in 1:(GridWidth*GridLength)) plot(test.som$codes[i,],type=l) Only when I use a larger grid (say for example 7*3 ) I get some of the representatives for the sin pattern. This must have something to do with the initialization of the grid, as the sin is so rare it is unlikely that I get it as a reference vector. Afterwards, because the selection for the training is also random it is also unlikely they are picked. I've been trying to modify some of the other parameters for the SOM also, but I would appreciatte some input to keep me going until I receive the reference books from my bookstore. Are my suspictions right? Should I be using the SOM for my study or should I look somewhere else? NOTE: I have no prior knowledge of whether the datasets I want to analyse will have rare cases or not or where they will be located. Thanks, Manuel __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] memory and swap space in ncdf
I've a linux system with 2Gb of memory which is not enough for reading a 446Mb netcdf file using ncdf: library(ncdf) ncold - open.ncdf(gridone.grd) Error: cannot allocate vector of size 1822753 Kb When I look at the free memory in my system I can see that none of the Swap space is being used by R. I am a newbie in linux and R, I've read the Memory help pages but still have some questions: can I use the swap space in R to solve my problem of lack of memory? if not, are there any ways to read the data apart from buying more RAM? Thanks, M __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] unlist kills R
When I try to unlist a very large list, R is killed without any other warning: A-as.list(as.data.frame(matrix(1:21639744,nrow=3578,ncol=6048))) with AA-unlist(A) or AA-c(A,recursive=TRUE) I get a R terminado (killed) and the end of the session. I think I'll need to get more RAM (now 1Gb, any other solutions welcomed) to be able to do this but, shouldn't I get a more gentle warning than the kill message? Thanks, Manuel platform i686-pc-linux-gnu arch i686 os linux-gnu system i686, linux-gnu status major2 minor0.1 year 2004 month11 day 15 language R __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] test multiple objects for being equal length
I could not find any help pages on How to test many objects for being of equal length Something like identical for more than two objects? x-1:6 y-1:10 z-3:5 ## For two objects I can do: identical(length(x),length(y)) ## For more than two I currently can do: length(unique(c(length(x),length(y),length(z==1 but there must be a better way. Thanks, M __ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] substitute accents
I have an openoffice spreadsheet with a column of character strings. Some of them contain accents. I want to read it in R so I have saved it as a csv file using Western Europe (ISO-8859-1) character set (the default, I've tried other sets but it doesn't help). R reads it fine with CharMatrix-read.csv(test.csv,header=F,sep=,,as.is=TRUE); Say I wan't to replace the 'o' with accent in the first cell. I've tried: gsub('ó','o', CharMatrix[1,1]) But, It doesn't make any substitution Trying to find a solution I input the character string in R and do the substitution: CharMatrix[1,1]-hóla gsub('ó','o', CharMatrix[1,1]) And it works. I think the difference is that when I now print the content of CharMatrix I get a \201 before the ó while I didn't get it with the openoffice imported csv file. I'm sure it is a problem with my understanding of how accents can be specified. Can someone give me any solutions / references? Thanks, M _ platform i686-pc-linux-gnu arch i686 os linux-gnu system i686, linux-gnu status major2 minor0.0 year 2004 month10 day 04 language R __ Renovamos el Correo Yahoo!: ¡100 MB GRATIS! Nuevos servicios, más seguridad __ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] substitute accents
$ locale LANG=en_GB LC_CTYPE=en_GB LC_NUMERIC=en_GB LC_TIME=en_GB LC_COLLATE=en_GB LC_MONETARY=en_GB LC_MESSAGES=en_GB LC_PAPER=en_GB LC_NAME=en_GB LC_ADDRESS=en_GB LC_TELEPHONE=en_GB LC_MEASUREMENT=en_GB LC_IDENTIFICATION=en_GB LC_ALL= $ locale charmap ISO-8859-1 I have tried changing the locales with no difference. Is this fine? And, no, I didn't get any warning message. My sistem is a debian sid under kde 3.3. Thanks, M --- Prof Brian Ripley [EMAIL PROTECTED] escribió: Can you please tell us what locale you are working in? This looks as if the problem might be the use of a UTF-8 locale, which R does not currently support and which some Linux distros have made their default. However, R does issue a warning -- so did you get one? On Thu, 25 Nov 2004, Manuel Gutierrez wrote: I have an openoffice spreadsheet with a column of character strings. Some of them contain accents. I want to read it in R so I have saved it as a csv file using Western Europe (ISO-8859-1) character set (the default, I've tried other sets but it doesn't help). R reads it fine with CharMatrix-read.csv(test.csv,header=F,sep=,,as.is=TRUE); Say I wan't to replace the 'o' with accent in the first cell. I've tried: gsub('ó','o', CharMatrix[1,1]) But, It doesn't make any substitution Trying to find a solution I input the character string in R and do the substitution: CharMatrix[1,1]-hóla gsub('ó','o', CharMatrix[1,1]) And it works. I think the difference is that when I now print the content of CharMatrix I get a \201 before the ó while I didn't get it with the openoffice imported csv file. I'm sure it is a problem with my understanding of how accents can be specified. Can someone give me any solutions / references? Thanks, M _ platform i686-pc-linux-gnu arch i686 os linux-gnu system i686, linux-gnu status major2 minor0.0 year 2004 month10 day 04 language R __ Renovamos el Correo Yahoo!: ¡100 MB GRATIS! Nuevos servicios, más seguridad __ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ Renovamos el Correo Yahoo!: ¡100 MB GRATIS! Nuevos servicios, más seguridad __ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html