Re: [R] Predictions with missing inputs
On Fri, 2011-02-11 at 20:51 -0500, Axel Urbiz wrote: Dear users, I'll appreciate your help with this (hopefully) simple problem. I have a model object which was fitted to inputs X1, X2, X3. Now, I'd like to use this object to make predictions on a new data set where only X1 and X2 are available (just use the estimated coefficients for these variables in making predictions and ignoring the coefficient on X3). Here's my attempt but, of course, didn't work. #fit some linear model to random data x=matrix(rnorm(100*3),100,3) y=sample(1:2,100,replace=TRUE) mydata - data.frame(y,x) mymodel - lm(y ~ ns(X1, df=3) + X2 + X3, data=mydata) summary(mymodel) #create new data with 1 missing input mynewdata - data.frame(matrix(rnorm(100*2),100,2)) mypred - predict(mymodel, mynewdata) Thanks in advance for your help! Axel. Axel, I think mice package solve your problem -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] can't find and install reshape2??
On Mon, 2010-10-04 at 10:27 +0930, Chris Howden wrote: Hi everyone, Im trying to install reshape2. But when I click on install package its not coming up!?!?! Im getting reshape, but no reshape2? Ive also tried download.packages(reshape2, destdir=c:\\) download.packages(Reshape2, destdir=c:\\)but no luck!!! Does anyone have any ideas what could be going on? Chris Howden Hi Chris, I have two guess: 1- You don't have installed 'stringr' pakage 2- Your R is outdated Try this two things and after this mail me -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] read columns of quoted numbers as factors
On Mon, 2010-10-04 at 09:39 -0700, james hirschorn wrote:
Suppose I have a data file (possibly with a huge number of columns), where
the
columns with factors are coded as 1, 2, 3, etc ... The default behavior
of
read.table is to convert these columns to integer vectors.
Is there a way to get read.table to recognize that columns of quoted numbers
represent factors (while unquoted numbers are interpreted as integers),
without
explicitly setting them with colClasses ?
__
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.
Hi James,
I think you solve ypur problem using the options colClasses in the
read.table command, something like this
rea.table('name.of.table',colClasses=c(rep(30,'integer'),rep(5,'numeric'),etc))
--
Bernardo Rangel Tura, M.D,MPH,Ph.D
National Institute of Cardiology
Brazil
__
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.
Re: [R] Calculating with tolerances
On Thu, 2010-09-09 at 09:16 +0430, Jan private wrote: Dear list, I am from an engineering background, accustomed to work with tolerances. For example, I have measured Q = 0.15 +- 0.01 m^3/s H = 10 +- 0.1 m and now I want to calculate P = 5 * Q * H and get a value with a tolerance +- What is the elegant way of doing this in R? Thank you, Jan Hi Jan, If I understood your problem this script solve your problem: q-0.15 + c(-.1,0,.1) h-10 + c(-.1,0,.1) 5*q*h [1] 2.475 7.500 12.625 -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] levene.test
On Wed, 2010-07-14 at 01:37 -0700, martanair wrote: I am trying to use Levene's test (of package car), but I do not understand quite well how to use it. '?levene.test' does not unfortunately provide any example. My data are in a data frame and correspond to 1 factor plus response. Could someone please give me an example about how to use the command levene.test(y, group) Thanks in advance, marta Hi Marta, levene.test is deprecated functions in car package as say the help of command (?levene.test). Instead this use leveneTest (?levene.test again)If you type ?leveneTest you get examples -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] F# vs. R
On Wed, 2010-07-07 at 17:31 +0200, Sergey Goriatchev wrote: Hello, everyone F# is now public. Compiled code should run faster than R. Anyone has opinion on F# vs. R? Just curious Best, S Sergey, F# is public, but is not open source. F# run in windows but run in AIX, linux, MAC, UNIX etc? Compiled code should run faster than R, but is precise? Compiled code should run faster than R, but is reliable? Compiled code should run faster than R, but have 2.440 packages for extend your capacities? Compiled code should run faster than R, but critical code is in C o FORTRAN So I think the F# is not a good alternative, if your concern is velocity dou you look Littler http://code.google.com/p/littler/ -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] How to determine if R is 64 bit compiled under Unix-alike?
On Mon, 2010-07-05 at 19:25 +0200, Przemek Grabowicz wrote: Under MacOS I had R64 executive and it was clear. Under Ubuntu, which I do not have administrative rights to, there is only R executive. It seems that I can allocate more than 3GB of memory, however not everything seems to work the same/right as with R64 under MacOS. Pms. Type .Machine$sizeof.pointer If respond is 8 your R is 64 bits -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Best way to compute a sum
On Sat, 2010-07-03 at 17:27 -0700, Roger Deangelis wrote: Hi Bernado, In many financial applications if you convert the dollars and cents to pennies( ie $1.10 to 110) and divide by 100 at the vary end you can get maintain higher precision. This applies primarily to sums. This is similar to keeping track of the decimal fractions which have exact representations in floating point. It is good to know where the errors come from. OK I suppose you could also improve the sum by understanding the decimal to binary rounding algorithms. I have noticed differences in computations between Sun hardware(FPUs) and Intel(FPUs). When presented the same set of operations, it appears that the floating point hardware do not do the operations in the same order. It's true Consider sprintf(%a, sum(rep(1.1,100))) # overestimates the sum 0x1.b8001p+6 0x1.b8001p+6 sprintf(%a, sum(rep(11,100))/10) # this gives the correct answer 0x1.b8p+6 110 = 10 1110 = 32 + 8 + 4 + 2 - 0x1.b8p+6 - tricky due to 53 bit and little endian (I think this is right) 0x1.b8p+6 note cmp - ifelse (sum(rep(1.1,100))==sum(rep(11,100))/10, equal, unequal) [1] unequal cmp - ifelse (sum(rep(1.1,100))sum(rep(11,100))/10, greater than, less than orequal) [1] greater than Hi Roger, First of all, is true sum(rep(1.1,100))==sum(rep(11,100))/10 [1] FALSE But is true too all.equal(sum(rep(1.1,100)),(sum(rep(11,100))/10)) [1] TRUE I think you need read about Guard Digits approach -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Best way to compute a sum
On Fri, 2010-07-02 at 21:23 -0700, Roger Deangelis wrote: Although it does not apply to your series and is impractical, it seems to me that the most accurate algorithm might be to add all the rational numbers whose sum and components can be represented without error in binary first, ie 2.5 + .5 or 1/16 + 1/16 + 1/8. You could also get very clever and investigate a sum that should have an exact binary representation when the individual components do not, ie .1 + .2 + .2 = .5 and correct the sum. Roger Roger I think you must read: What Every Computer Scientist Should Know About Floating-Point Arithmetic ( http://docs.sun.com/source/806-3568/ncg_goldberg.html ) I think your question and others like this question is answer in this paper -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] identifying odd or even number
On Thu, 2010-07-01 at 08:40 -0700, Yemi Oyeyemi wrote: Hi, I run into problem when writing a syntax, I don't know syntax that will return true or false if an integer is odd or even. Thanks OYEYEMI, Gafar Matanmi Department of Statistics University of Ilorin Hi Yemi! Your problem is solve wiht ! and %%, Look this example x - 1:11 x [1] 1 2 3 4 5 6 7 8 9 10 11 x%%2 [1] 1 0 1 0 1 0 1 0 1 0 1 !(x%%2) [1] FALSE TRUE FALSE TRUE FALSE TRUE FALSE TRUE FALSE TRUE FALSE -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Calculation of r squared from a linear regression
On Fri, 2010-06-11 at 01:16 -0700, Sandra Hawthorne wrote: Hi, I'm trying to verify the calculation of coefficient of determination (r squared) for linear regression. I've done the calculation manually with a simple test case and using the definition of r squared outlined in summary(lm) help. There seems to be a discrepancy between the what R produced and the manual calculation. Does anyone know why this is so? What does the multiple r squared reported in summary(lm) represent? # The test case: x - c(1,2,3,4) y - c(1.6,4.4,5.5,8.3) dummy - data.frame(x, y) fm1 - lm(y ~ x-1, data = dummy) summary(fm1) betax - fm1$coeff[x] * sd(x) / sd(y) # cd is coefficient of determination cd - betax * cor(y, x) Thanks. Sorry Sandra, But the problem in yours script. Look this x - c(1,2,3,4) y - c(1.6,4.4,5.5,8.3) dummy - data.frame(x, y) fm1 - lm(y ~ x, data = dummy) summary(fm1) Call: lm(formula = y ~ x, data = dummy) Residuals: 1 2 3 4 -0.17 0.51 -0.51 0.17 Coefficients: Estimate Std. Error t value Pr(|t|) (Intercept) -0.3500 0.6584 -0.532 0.6481 x 2.1200 0.2404 8.818 0.0126 * --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 0.5376 on 2 degrees of freedom Multiple R-squared: 0.9749, Adjusted R-squared: 0.9624 F-statistic: 77.76 on 1 and 2 DF, p-value: 0.01262 betax - fm1$coeff[2] * sd(x) / sd(y) # cd is coefficient of determination cd - betax * cor(y, x) cd x 0.974924 The formula fm1$coeff[2] * sd(x) / sd(y) is valid only the model have a intercept... -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] How sample without replacement on more than one variables?
On Sun, 2010-05-23 at 00:56 -0700, dusadrian wrote: This might help, depending on your exact needs: v1 - sample(letters[1:2], 10, replace=TRUE) v2 - sample(letters[3:4], 10, replace=TRUE) v3 - sample(letters[5:6], 10, replace=TRUE) aa - data.frame(v1=v1, v2=v2, v3=v3) And now is simple, sample the row of data frame aa[sample(1:nrows(aa),3),] -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Minimization problem
On Thu, 2010-05-20 at 01:35 -0700, Fred wrote:
Dear R users,
I am trying to minimize two function simultaneously in R,
function(x)
minimize x[1],x[2],x[3]
mean(distribution(x1,x2,x3) ) - observed mean
std(distribution(x1,x2,x3)) - observed std
What I want to achieve is that simulated mean and standard deviation
of distribution related to x1 x2 x3 would be close to observed mean
and observed standard deviation.
is there any function in R can reach this?
Thank you for the help first .
F.
Hi!
Do you need use optim, something like this
test - function(parameters){
m.error - mean(distribution(x1,x2,x3) ) - observed mean
m.sd - std(distribution(x1,x2,x3)) - observed std
res - cbind(m.error,sd.error)
return(res)
}
--
Bernardo Rangel Tura, M.D,MPH,Ph.D
National Institute of Cardiology
Brazil
__
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.
Re: [R] Simple question on binning data
On Thu, 2010-05-13 at 17:31 -0400, Carl Witthoft wrote: It's very simple to write a binit() function. If all you want to do is e.g., bin 107 values into sums of 10 at a time, then write a loop that sums x[10*i:11*i-1] (not tested and not syntactically correct). The one I wrote for myself discards any partial bin (101-107 in my example) and leaves a warning note that this took place. Carl Hi Carl, I think the syntactically correct is x[10*i:(11*i-1)] -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Generate Numbers
On Sun, 2010-04-04 at 04:03 -0700, kende jan wrote: Hello, How can I generate randomly in R a sample of skewed data with first quartile is 540 and third quartile is 715. I need a sample of 100 cases. Thank you for your help Jan Hi Jan, Sorry my direct approach but why you don't use data - c(runif(25,number,540),runif(50,540,715),runif(25,715,number)) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] compare two fingerprint images
On Fri, 2010-04-02 at 20:52 +0200, Juan Antonio Gil Pascual wrote: Hello I wanted to compare two fingerprint images. How do you do with R?. Is there a role for cross-correlation of images? Thanks Hi juan, You was a quite vage in your question, so I don't know exactly what you need but try this; require(ReadImages) x - read.jpeg(image1) x1 - read.jpeg(image2) table(x1==x) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] compare two fingerprint images
On Sat, 2010-04-03 at 21:18 +0200, Juan Antonio Gil Pascual wrote: Hi Bernado I need to compare two fingerprint images and let me know if you can do with R. I have used the technique of minutiae but it seems to work better with the cross-correlation and wanted to know if you can do with R. Thank you very much Juan Bernardo Rangel Tura escribió: On Fri, 2010-04-02 at 20:52 +0200, Juan Antonio Gil Pascual wrote: Hello I wanted to compare two fingerprint images. How do you do with R?. Is there a role for cross-correlation of images? Thanks Hi juan, You was a quite vage in your question, so I don't know exactly what you need but try this; require(ReadImages) x - read.jpeg(image1) x1 - read.jpeg(image2) table(x1==x) Juan, I don't know cross-correlation to images but I know this for time series in this case: ccf(ts(as.numeric(x)),ts(as.numeric(x1))) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Changing content of column in data.frame + efficient join extraction between 2 data.frames
On Tue, 2010-03-23 at 09:37 +, Corrado wrote: Dear R users, I have 2 SpatialPointsDataFrame's, pcs and East. The column str_1 in the first (pcs) is: pcs[0:4,] coordinates cat str_1 int_1 int_2dbl_1 dbl_2 1 (101000, 263000) 1 SM06B 101000 263000 4.978915 -4.293668 2 (101000, 265000) 2 SM06C 101000 265000 4.960478 -4.266742 3 (101000, 267000) 3 SM06D 101000 267000 4.912984 -4.246849 4 (101000, 269000) 4 SM06E 101000 269000 4.613309 -4.185405 The column str_1 in the second (East) is: East[0:4,] coordinates str_1 1 (489000, 215000) sp81x 2 (489000, 217000) sp81y 3 (493000, 209000) sp90j 4 (495000, 209000) sp90p I would like to do 2 things: 1) I would like to change the format of the column str_1 in the first to be the same that it is in the second, that is I need to remove the inverted commas and I need to make it lower case. 2) I would like to extract the rows from the first one (pcs) where pcs$str_1 is the same as East$str_1. I have even tried regexp, but cannot modify the content of pcs$str_1 to remove the inveretd commas and change the case to lowercase. How do I do that? Regards Hi Corrado! First: tolower(pcs$str_1) change to lower case Second: try merge (East,pcs,by.x=str_1,by.y=str_1) to fusion data frames Third: I don't recreate your database East in my computer do you give a small part to I try solve your problem? -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] queue simulation
On Mon, 2010-03-22 at 02:49 -0600, Carlos Ernesto Lopez Nataren wrote: Hello everybody :) I am trying to simulate a queue with times of arrival to the queue and time taken to dispatch every member of the queue coming from two exponential distributions, I am interested in knowing the number of people at any time and the time that takes every member of this queue to be dispatched. I thought this was gonna be an easy task but I've failed to try to simulate this, is there any package that does this already? Any help will be greatly appreciated. Thank you very much Carlos Carlos, if I understand your problem do you need know the time for each person in a queue dispatched. You say this time fit a exponential distribution, so do you have a rate of dispatch (x dispatch per time) If I want generate 10 times of dispatch rate 0.4, use the command times - rexp(10,0.4) If I need the total delay for each person, use the command cumsum(times) If I need the average time in the queue, use the command means(cumsum(times)) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] different forms of nls recommendations
On Sat, 2010-03-20 at 14:55 -0800, emorway wrote:
Hello,
Using this data:
http://n4.nabble.com/file/n1676330/US_Final_Values.txt US_Final_Values.txt
and the following code i got the image at the end of this message:
US.final.values-read.table(c:/tmp/US_Final_Values.txt,header=T,sep= )
US.nls.1-nls(US.final.values$ECe~a*US.final.values$WTD^b+c,data=US.final.values,start=list(a=2.75,b=-0.95,c=0.731),trace=TRUE)
f.US1-function(x){coef(US.nls.1)[a]*x^coef(US.nls.1)[b]+coef(US.nls.1)[c]}
xvals.US1-seq(min(US.final.values$WTD),max(US.final.values$WTD),length.out=75)
yvals.US1-f.US1(xvals.US1)
Rsq.nls.1-sum((predict(US.nls.1)-mean(US.final.values$ECe))^2/sum((US.final.values$ECe-mean(US.final.values$ECe))^2))
plot(US.final.values$WTD,US.final.values$ECe,col=red,pch=19,cex=.75)
lines(xvals.US1,yvals.US1,col=blue)
but the r^2 wasn't so hot.
Rsq.nls.1
[1] 0.2377306
So I wanted to try a different equation of the general form a/(b+c*x^d)
US.nls.2-nls(US.final.values$ECe~(a/(b+c*US.final.values$WTD^d)),data=US.final.values,start=list(a=100.81,b=73.7299,c=0.0565,d=-6.043),trace=TRUE,algorithm=port)
but that ended with Convergence failure: false convergence (8). I tried
Hi emorway,
Do you have 657 obs and 4 parameters to fit.
In my opinion you have few obs...
I think do you fit in steps:
US.nls.2-nls(ECe~(a/(b + c *
WTD^d)),data=US.final.values,start=list(a=100.81,b=73.7299,c=0.0565,d=-6.043),trace=TRUE,algorithm=port)
temp_nls1 - nls(ECe~(100/(73 + .05 *
WTD^d)),data=US.final.values,start=list(d=-6.043),trace=TRUE,algorithm=port)
temp_nls2 - nls(ECe~(100/(73 + .05 *
WTD^d)),data=US.final.values,start=list(d=-1.01613),trace=TRUE,algorithm=port)
temp_nls3 - nls(ECe~(100/(73 + c *
WTD^(-1.01613))),data=US.final.values,start=list(c=0.05),trace=TRUE,algorithm=port)
temp_nls4 - nls(ECe~(100/(73 + c *
WTD^(-1.01613))),data=US.final.values,start=list(c=-14.7127),trace=TRUE,algorithm=port)
temp_nls5 - nls(ECe~(100/(b-14.7127 *
WTD^(-1.01613))),data=US.final.values,start=list(b=73),trace=TRUE,algorithm=port)
temp_nls6 - nls(ECe~(100/(b-14.7127 *
WTD^(-1.01613))),data=US.final.values,start=list(b=70.4936),trace=TRUE,algorithm=port)
temp_nls7 - nls(ECe~(a/(70.4936-14.7127 *
WTD^(-1.01613))),data=US.final.values,start=list(a=100),trace=TRUE,algorithm=port)
0: 2243.9898: 100.000
1: 2122.8218: 106.219
2: 1359.8819: 187.530
3: 1359.8819: 187.530
--
Bernardo Rangel Tura, M.D,MPH,Ph.D
National Institute of Cardiology
Brazil
__
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.
Re: [R] fisher.test gives p1
On Thu, 2010-03-04 at 11:15 -0500, Jacob Wegelin wrote: The purpose of this email is to (1) report an example where fisher.test returns p 1 (2) ask if there is a reliable way to avoid p1 with fisher.test. If one has designed one's code to return an error when it finds a nonsensical probability, of course a value of p1 can cause havoc. Example: junk-data.frame(score=c(rep(0,14), rep(1,29), rep(2, 16))) junk-rbind(junk, junk) junk$group-c(rep(DOG, nrow(junk)/2), rep(kitty, nrow(junk)/2)) table(junk$score, junk$group) DOG kitty 0 1414 1 2929 2 1616 dput(fisher.test(junk$score, junk$group)$p.value) 1.12 Hi jacob, I think this is cover in R FAQ 7.31, but look this command all.equal(dput(fisher.test(matrix(c(14,14,29,29,16,16),byrow=T,ncol=2))$p.value),1) 1.12 [1] TRUE P.S R FAQ 7.31 - http://cran.r-project.org/doc/FAQ/R-FAQ.html#Why-doesn_0027t-R-think-these-numbers-are-equal_003f -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Problem in performing goodness of fit test in R.
On Sun, 2010-02-14 at 12:42 +0500, Faiz Rasool wrote: I am trying to perform goodness of fit test using R. I am using this website http://wiener.math.csi.cuny.edu/Statistics/R/simpleR/stat013.html for help. However, I am unable to carry out the test successfully. My code follows. It is taken from the website just mentioned. freq=c(22,21,22,27,22,36) # frequencies obtained after rolling the dice 150 times. prob=c(1,1,1,1,1,1)/6 # specify expected frequency for each category. chisq.test(freq,p=prob) # I do not know what this line means. I just followed instructions on the website. The erorr I receive is erorr in chisq.test(freq,p=prob)/6 probabilities must sum to 1 I am very new to R, so any help would be appreciated. Faiz. Faiz, Well ... In my computer( Phenom X4 9650, runing Ubuntu 9.10 and R 2.10.1) the script work freq=c(22,21,22,27,22,36) # frequencies obtained after rolling the dice 150 times. prob=c(1,1,1,1,1,1)/6 # specify expected frequency for each category. chisq.test(freq,p=prob) # I do not know what this line means Chi-squared test for given probabilities data: freq X-squared = 6.72, df = 5, p-value = 0.2423 About the third line You must read ?chisq.test for better know the command, but you execute one chi-square test with uniform probability distribution -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] using tabble
On Mon, 2010-01-18 at 17:17 +0800, Fan Dan wrote: Dear Helper: I am using R to finish a project which is near the dealine. Please help me ! I need to construct a table. The name of the table is the sequence from 1 to 7 with some missing values arbitrarily. For example 1 to 7 12 36 7 12 3 45 5 7 I want to add the missing numbers of the title which should be seq(1,7) In this case, 4 and 5 and the corresponding data is 0 So after revised the table should be like 12 3 4 5 6 7 12 3 45 0 0 5 7 Please tell me how I really appreaciate your time and helpn Yours Wolfgang AMadeurs [[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. Hi! If i understood your problem is easy create a fake data for example: set.seed(123) x-rep(sample(1:7,5),3) table(x) x 3 4 5 6 7 3 3 3 3 3 In this example the values 1 and 2 are missing Solution change your variable to factor and determine your leveal x.factor-factor(x,levels=1:7) table(x.factor) x 1 2 3 4 5 6 7 0 0 3 3 3 3 3 -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Help with function fitdistr in MASS
On Sat, 2010-01-02 at 23:20 -0800, Saji Ren wrote: Hi, R users: I want to fit my data into a normal distribution by using the command fitdistr in MASS. I changed my data class from ts to numeric by class(mydata)=numeric but after using fitdistr, I got the result below fitdistr(mydata,normal) meansd NA NA (NA) (NA) the help doc of fitdistr does not mention anything about that, thus I need your help. Thank you in advanced, Saji from Shanghai Hi Sajj, You hava NA in your data try: fitdistr(na.exclude(mydata),normal) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Problem with Cannot compute correct p-values with ties
On Wed, 2009-12-02 at 16:52 +0800, Zhijiang Wang wrote: Dear All, 1. why did the problem happen? 2. How to solve it? -- Best wishes, Zhijiang Wang Well... The algorithm for Mann-whitney test have problem with ties To solve you can use jitter a-1:10 b-1:10 wilcox.test(a,b) Wilcoxon rank sum test with continuity correction data: a and b W = 50, p-value = 1 alternative hypothesis: true location shift is not equal to 0 Warning message: In wilcox.test.default(a, b) : cannot compute exact p-value with ties wilcox.test(a,jitter(b)) Wilcoxon rank sum test data: a and jitter(b) W = 49, p-value = 0.9705 alternative hypothesis: true location shift is not equal to 0 look ?jitter for more information -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] areg (stata) equivalent in R?
On Mon, 2009-11-30 at 22:20 -0500, John K. Williams wrote: hi, i am looking to reproduce a study done in stata in R, where a regression was done while absorbing a categorical variable. i am new to R, i've i installed the design package but haven't been able to find an applicable function. thanks for any help. John, Do you show a example for this command? I don't know stata so I don't help you -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] subset or condition as argument to a function
On Fri, 2009-11-20 at 17:40 -0800, Santosh wrote:
Dear Rxperts!
I was wondering if it is possible to write a function which can take in
argument of a subset or condition.. Of course, I am aware of the alternate
methods like coplot, par.plot, xyplot etc... I am specifically interested in
using conditions/subsets with plot..
A simple fragmented example is shown here..
pltit - function(y,x,dat,dat1,dat2,sbst) {
plot(y~x, data=dat, subset=sbst)
lines(y~x,data=dat1, subset=subst)
points(y~x,data=dat2,subset=subst)
}
pltit(profit,weeks,dat=zone1, sbst='group==1sales200')
Could you also suggest pointers/references/examples on efficient ways to
plot simulated data overlaid with observed data?
Have a good weekend!
Hi Santosh,
IMHO your script is
pltit - function(y,x,dat,dat1,dat2,sbst) {
subs - subset(dat,sbst)
with(subs,plot(y~x))
subs1 - subset(dat1,sbst)
with(subs1,lines(y~x))
subs2 - subset(dat2,sbst)
with(subs2,points(y~x))
}
--
Bernardo Rangel Tura, M.D,MPH,Ph.D
National Institute of Cardiology
Brazil
__
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.
Re: [R] negative log likelihood
On Sun, 2009-11-08 at 11:12 -0800, mat7770 wrote: I have two related variables, each with 16 points (x and Y). I am given variance and the y-intercept. I know how to create a regression line and find the residuals, but here is my problem. I have to make a loop that uses the seq() function, so that it changes the slope value of the y=mx + B equation ranging from 0-5 in increments of 0.01. The loop also needs to calculate the negative log likelihood at each slope value and determine the lowest one. I know that R can compute the best regression line by using lm(y~x), but I need to see if that value matches the loop functions. Hi If I understand your question you need extract log-likelihood for a linear model, so you need using logLik command, see example: set.seed(1) x-rpois(16,6) y-2*x+3+rnorm(16,sd=3) model-lm(y~x) logLik(model) 'log Lik.' -40.1177 (df=3) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] operation with if/else on a dataframe
On Thu, 2009-10-29 at 01:47 -0700, Fran100681 wrote: Hi to all, I have this dataframe (I show the first six rows) head(table) A RFold.Change P.Value Count1 Count2 1 ENSRNOE002_at 0 1.13 0.601 1 2 ENSRNOE009_at 0-1.04 0.733 3 3 ENSRNOE020_at 0-1.08 0.680 0 4 ENSRNOE021_at 0-1.31 0.201 2 5 ENSRNOE023_at 0-1.06 0.643 3 6 ENSRNOE024_at 0-1.14 0.403 3 I would like to generate a function that determine for each row a new value (resulting in a new vector of values to add to the dataframe). The function should give for every row the same value showed in column R. However,I need of this because not all the R-values reported in this table are correctly determined following the criteria mentioned below. These new values calculated by the function must be: 1)UP if all the following conditions are verified in a certain row: Fold.Change is = +1.5, P.Value is 0.05 Count1 = 2 2)DOWN if all the following conditions are verified in a certain row: Fold.Change is = -1.5, P.Value is 0.05 Count2 = 2 3) 0 if both of previous conditions are not verified So I have set these fllowing parameters (new objects) because I'll have to repeat this procedure to different dataframes in which the order of columns of interest might change (So I can change these parameters depending on the order of the columns in any different table) Fold.change - 3 #(because in this table, Fold.Change value is the third column and so on...) P.Value - 4 Count1 - 5 Count2 - 6 [ Quote text] This is my problem, I cannot use the function to recalculate values in R-column for all rows in my dataframe. I don't understand where is the problem, can someone help me? Thanks a lot!! Francesco Francesco, I think you solve this problem with a simple way. Remember in R the most function and operations are vectorized so look this example: set.seed(123) x-rpois(20,5) y-rpois(20,15) z-rpois(20,10) dta-data.frame(x,y,z) dta dta$NEW-ifelse(x5 y15 z10,UP, ifelse(x5 y15 z10,DOWN, 0)) dta First, I use ifelse command to simplify your nested conditional situation. Second, I know that R test this nested condition in order so the first position will result test x[1],y[1] and z[1], the second postion will result test x[2],y[2] and z[2] ... The new vector result is the same order the original data.frame so I use dta$NEW to create a new column in data.frame -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] attach
On Wed, 2009-10-14 at 07:21 +0200, Christophe Dutang wrote: Hi all, I have a question regarding the memory usage for the attach function. Say I have a data.frame inputdat that I create with read.csv. I would like to know what happens on the memory side when I use attach(inputdata) Is there a second allocation of memory for inputdata? Then I'm using eval on a expression which depends on the columns of inputdata. Is it better not to use attach function? Thanks in advance Christophe Well, if you attach a data.frame twice times, it use your memory twice times. I don't use attach I prefer with(data.frame, command) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] stochastic process
On Tue, 2009-10-13 at 14:43 +0800, 刘哲 wrote: Hi, I'm a student in China, and I never used R before. I'm now wondering how to simulate a sample of Markov chain of ,say 500 entries with a certain transition matrix. Thanks a lot. Hi Well the naive script for you problem is something like this # # initial condition # initial-matrix(c(300,200,0),nrow=1) initial [,1] [,2] [,3] [1,] 300 2000 # # transition matrix (third state is absortive) # transition-matrix(c( .5,.3,.2, .3,.4,.3, 0,0,1 ),nrow=3,byrow=T) transition [,1] [,2] [,3] [1,] 0.5 0.3 0.2 [2,] 0.3 0.4 0.3 [3,] 0.0 0.0 1.0 # # first step # S1-initial%*%transition S1 [,1] [,2] [,3] [1,] 210 170 120 # # Second step # S2-S1%*%transition S2 [,1] [,2] [,3] [1,] 156 131 213 If you use a large number of step is more pratical use command mtx.exp of Biodem package, something like this require(Biodem) # # 10th step direct # initial%*%mtx.exp(transition,10) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] row selection
On Thu, 2009-10-08 at 16:14 -0400, Ashta wrote:
Hi all,
I have a matrix named x with N by C
I want to select every 5 th rrow from matrix x
I used the following code
n- nrow(x)
for(i in 1: n){
+ b - a[i+5,]
b
}
Error: subscript out of bounds
Can any body point out the problem?
Hi Ashta,
If I understand your request you need select row 5,10,15, ...
In this case you can use this script:
x[1:nrow(n)%%5==0]
--
Bernardo Rangel Tura, M.D,MPH,Ph.D
National Institute of Cardiology
Brazil
__
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.
Re: [R] Setting a mirror permanently on R on ubuntu
On Sat, 2009-10-10 at 04:54 +0100, Lazarus Mramba wrote: Dear all, I seem to have many problems as I run R on my ubuntu system. want to set a mirror so that anytime I use the command install.packages, it does not ask me for which mirror to use but go direct. This is because of the error I keep on getting below and I dont know how to solve it. Please help. Kind regards, Lazarus Lazarus, use options(repos=mirror URL) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Satellite ocean color palette?
On Fri, 2009-10-09 at 11:51 -0700, Tim Clark wrote: Dear List, Is there a color palette avaliable similar to what is used in satellite ocean color imagery? I.e. a gradient with blue on one end and red on the other, with yellow in the middle? I have tried topo.colors(n) but that comes out more yellow on the end. I am looking for something similar to what is found on the CoastWatch web page: http://oceanwatch.pifsc.noaa.gov/imagery/GA2009281_2009282_sst_2D_eddy.jpg Thanks! Tim Tim Clark Department of Zoology University of Hawaii Tim, You can make a palette in R, using colorRampPalette, look this example Satelite.Pallete - colorRampPalette(c(blue3,cyan,aquamarine,yellow,orange,red)) require(fields) image.plot(volcano, col = Satelite.Pallete(500), legend.lab=Scale) contour(volcano, levels = seq(90, 200, by = 5), add = TRUE) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] How to use Subpopulation data?
On Fri, 2009-10-02 at 00:42 +, KABELI MEFANE wrote: Thanks But it seems like you don't get my problem.Do you mean that there is something wrong with the code as it seems like what you are doing is suggesting different ways to write a code. Will i get to use the variable that have been name previously like if i want to calculate the standard deviation of stratum hypermarket in a sample not population, the first start is to check if it would help me is to check the length() of different levels if it is not same as that of the original data Best Regards I rest my case, i might dream it Sorry I don't understand your problem early. In your case you have 2 data.frames: sampleframe and sample.strat so you need show to R where calculate sd or length or summary set.seed(1)# only to make this example don't use in your code stra-strata(sampleframe,c(type),size=c(20,80,200,300,400),method=srswor) sample.strat-getdata(sampleframe,stra) length of value in sampleframe where type is Hypermarket: length(sampleframet$value[sampleframe$type==Hypermarket]) [1] 100 SD of value in sampleframe where type is Hypermarket: sd(sampleframe$value[sampleframet$type==Hypermarket]) [1] 4586.854 length of value in sample.strat where type is Hypermarket: length(sample.strat$value[sample.strat$type==Hypermarket] ) [1] 20 SD of value in sample.strat where type is Hypermarket: sd(sample.strat$value[sample.strat$type==Hypermarket] ) [1] 4679.336 -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] aproximate a titration kurve to the measure data.
On Thu, 2009-10-01 at 10:18 -0700, awayguy wrote: yes, halo thank you. my measure data: v2 - c(0, 2, 4, 6, 6.2, 6.4, 6.6, 6.8, 7, 7.2, 7.4, 7.6, 7.8, 8, 10, 12, 14) ph2 - c(12.10, 11.94, 11.68, 11.11, 10.91, 10.74, 10.47, 9.71, 7.1, 4.24, 3.3, 3.08, 2.98, 2.86, 2.33, 2.11, 1.98) with regards Ok temp-data.frame(ph2,v2) drm(temp,, fct = LL.4()) A 'drc' model. Call: drm(formula = temp, fct = LL.4()) Coefficients: b:(Intercept) c:(Intercept) d:(Intercept) e:(Intercept) 39.420 2.425 11.487 7.002 -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] aproximate a titration kurve to the measure data.
On Wed, 2009-09-30 at 10:55 -0700, awayguy wrote: Halo i'm studying chemistry, today we made an experiment and i have to draw a titration kurve for my mess data. we can do it on a mm paper, or we can also use a programe. people from chemistry recomend R last year i studied civil eng. and we used Matlab, as I see, R ist very similar to it, but its got other comands. But i think R would be a good help for some exercises. so my main question is: i have some measurement data from my titration, and I want aproximate a kurve to this data. is it possible to do it with R? a titration kurve looks like this: http://www.nabble.com/file/p25685986/acetic-acid-titration-curve.png hope you can help me, and yes when its possible, if you know something like a tutorial then i would be glad if you could post it. with regards Hi, I think do you need use drc package: drc: Analysis of dose-response curves Analysis of one or multiple curves with focus on concentration-response, dose-response and time-response curves used, for example in biology, environmental sciences, medicine, pharmacology, toxicology. -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] How to use Subpopulation data?
On Thu, 2009-10-01 at 10:06 +, KABELI MEFANE wrote: Dear Helpers I have a sample frame and i have sampled from it using three methods and now i want to calculate the statistics but i only get the population parameters. H - matrix(rnorm(100, mean=5, sd=5000)) sampleframe=data.frame(type=c(rep(H,100)),value=c(H)) sampleframe str=strata(sampleframe,c(type),size=c(20,), method=srswor) Try using str=strata(sampleframe,c(type),size=20, method=srswor) or better str - strata(sampleframe,c(type),size=20, method=srswor) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] How to use Subpopulation data?
On Thu, 2009-10-01 at 13:34 +, KABELI MEFANE wrote: ## package sampling stra=strata(sampleframe,c(type,value,rating),size=c(20,80,200,300,400), method=srswor) sample.strat-getdata(sampleframe,stra) sample.strat Try: stra-strata(sampleframe,size=c(20,80,200,300,400),method=srswor) sample.strat-getdata(sampleframe,stra) sample.strat -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] PCA or CA
On Tue, 2009-09-29 at 17:02 +, Paul Dennis wrote: Dear all I have a data set for which PCA based between group analysis (BGA) gives significant results but CA-BGA does not. I am having difficulty finding a reliable method for deciding which ordination technique is most appropriate. I have been told to do a 1 table CA and if the 1st axis is2 units go for CA if not then PCA. Another approach is that described in the Canoco manual - perform DCA and then look at the length of the axes. I used decorana in vegan and it gives axis lengths. I assume that these are measured in SD units. Anyway the manual say if the axis length is 3 go for PCA,4 use CA and if intermediate use either. Are either of these approaches good/valid/recommended or is there a better method? Thanks Paul Hi Paul I think that Ca is Correspondence Analysis and PCA is Principal Component Analysis, right? In this case, if all variables is numeric do you must use PCA, if all variables is factor do you must use CA. If you have a mixed variables do you have a problem, in most case I convert numeric variables to factor (with loss of information) and make CA -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Color index in image function
On Thu, 2009-09-10 at 05:27 -0700, FMH wrote: Thank you for the scripts, but the label and the values in the x and y-axis suddently dissapear even the 'axis' function is used as stated in the command below. Could you help on this? axis(1, at = seq(100, 800, by = 100)) axis(2, at = seq(100, 600, by = 100)) How could i add a tittle on top of color index? Thank you Fir Try this Brazilan.Pallete - colorRampPalette(c(green,yellow, blue)) require(fields) image.plot(volcano, col = Brazilan.Pallete(50), legend.lab=Scale) contour(volcano, levels = seq(90, 200, by = 5), add = TRUE) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Color index in image function
On Wed, 2009-09-09 at 02:33 -0700, FMH wrote: Thank you for the hints, but how could i add the grid lines which have numbers, representing the height of the volcano on the image. Thank you So I think this script is what you need Brazilan.Pallete - colorRampPalette(c(green,yellow, blue)) require(fileds) image.plot(volcano, col = Brazilan.Pallete(50), axes = FALSE) contour(volcano, levels = seq(90, 200, by = 5), add = TRUE) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil - Original Message From: Bernardo Rangel Tura t...@centroin.com.br To: FMH kagba2...@yahoo.com Sent: Tuesday, September 8, 2009 10:14:07 AM Subject: Re: [R] Color index in image function On Mon, 2009-09-07 at 07:59 -0700, FMH wrote: Thank you for the tips. I have manage to run your script, but was still never get the way to include the color index beside the image which could explain the intensity of the color from the lower index(green) to the higher index(blue). This color index might be represented by an increasing of color index in another table beside the image, started from green followed by green-yellow, yellow, yellow-blue and blue? Could someone please advice on this matter? Cheers Fir Hi FHM, Well If you desire one color index in a imageplot I don't know solve your problem. But in your scirptyou use image and image(x, y, volcano, col = terrain.colors(100), axes = FALSE) contour(x, y, volcano, levels = seq(90, 200, by = 5), add = TRUE, col = peru) In this case I suggest you use Brazilan.Pallete - colorRampPalette(c(green,yellow, blue)) filled.contour(volcano, color = Brazilan.Pallete) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Color index in image function
On Sat, 2009-09-05 at 04:14 -0700, FMH wrote: Dear All, I was looking for the color index in image function, such as from topo.colors(n) and etc. but still never found it. For instance, from the help menu. ### # Volcano data visualized as matrix. Need to transpose and flip # matrix horizontally. image(t(volcano)[ncol(volcano):1,]) # A prettier display of the volcano x - 10*(1:nrow(volcano)) y - 10*(1:ncol(volcano)) image(x, y, volcano, col = terrain.colors(100), axes = FALSE) contour(x, y, volcano, levels = seq(90, 200, by = 5), add = TRUE, col = peru) axis(1, at = seq(100, 800, by = 100)) axis(2, at = seq(100, 600, by = 100)) box() title(main = Maunga Whau Volcano, font.main = 4) # From the script above, it yields a beautiful image of volcano with variety of colors but i have to list down the color index that could show the meaning of each color in my thesis. Could someone please help me to extract this color index? Thank you Fir If I understand your question you need change the Palette of image plot. So you need use colorRampPalette look my example Brazilan.Pallete - colorRampPalette(c(green,yellow, blue)) image(x, y, volcano, col = Brazilan.Pallete(50), axes = FALSE) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] correlation between two 2D point patterns?
On Sun, 2009-08-30 at 07:51 +0100, William Simpson wrote: Suppose I have two sets of (x,y) points like this: x1-runif(n=10) y1-runif(n=10) A-cbind(x1,y1) x2-runif(n=10) y2-runif(n=10) B-cbind(x2,y2) I would like to measure how similar the two sets of points are. Something like a correlation coefficient, where 0 means the two patterns are unrelated, and 1 means they are identical. And in addition I'd like to be able to assign a p-value to the calculated statistic. cor(x1,x2) cor(y1,y2) gives two numbers instead of one. cor(A,B) gives a correlation matrix I have looked a little at spatial statistics. I have seen methods that, for each point, search in some neighbourhood around it and then compute the correlation as a function of search radius. That is not what I am looking for. I would like a single number that summarises the strength of the relationship between the two patterns. I will do procrustes on the two point sets first, so that if A is just a rotated, translated, scaled, reflected version of B the two patterns will superimpose and the statistic I'm looking for will say there is perfect correspondence. Thanks very much for any help in finding such a statistic and calculating it using R. Bill Hi Bill, If your 2 points set is similar I expect your Euclidean distance is 0, so I suggest this script: dist-sqrt((x1-x2)^2+(y1-y2)^2) # Euclidean distance t.test(dist) # test for mean equal 0 -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] SQL on R
On Sat, 2009-07-04 at 22:24 -0700, JoK LoQ wrote: I'm dealing with lots of columns and conditions, wats t best way to deal with that? How do I work with SQL on R? the manual is quite confuse talking about that. Do I need a package? I don't understand your question, but if you think use SQL in your data frames do you use SQLDF (http://cran.r-project.org/web/packages/sqldf/index.html) if you thinks use a database server and access it in R i sugest you use RMySQL (http://cran.r-project.org/web/packages/RMySQL/index.html) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Convert a lis to matrix
On Sun, 2009-06-07 at 00:14 -0400, Manisha Brahmachary wrote: Hello, This is an urgent request. I want to convert a list of 3 elements into a matrix and I am not sure how to do it. The list looks like this: List of 3 $ : num [1:15364, 1] 0.133 0.622 0.588 1.024 0.583 ... ..- attr(*, dimnames)=List of 2 .. ..$ : chr [1:15364] 6420681 3610072 2260458 60689 ... .. ..$ : NULL $ : num [1:15364, 1] 3.159 0.265 0.522 1.912 3.380 ... ..- attr(*, dimnames)=List of 2 .. ..$ : chr [1:15364] 6420681 3610072 2260458 60689 ... .. ..$ : NULL $ : num [1:15364, 1] 3.214 0.277 1.447 1.827 2.054 ... ..- attr(*, dimnames)=List of 2 .. ..$ : chr [1:15364] 6420681 3610072 2260458 60689 ... .. ..$ : NULL I'm not sure if understood your question, but look this code a-list(B=1:3,C=2:4,D=3:5) matrix(unlist(a),ncol=3) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
[R] Help with a cumullative Hazrd Ratio plot
Hi R-masters I need help to make modified cumulative hazard ratio plot. I need create a common plot but with the number of subjects in risk each ticks times for two different groups in bottom of plot (I put one example in attach). Do you know a routine for this? Is possible create a routine for this? In this case with how commands? Thanks in advance! -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] R and McAfee 8.5
On Mon, 2009-05-11 at 09:54 -0300, Joyce, Warren wrote: Hi, I have been working with R for the last year and using the UKFSST package to look at satellite tag track data and SST information. Fpr those not familiar, the package uses the positions estimated by the satellite tags themselves and the associated SST data from servers (in this case, from the University of Hawaii and the NOAA Coastwatch website) for the time preiod to estimate a reasonable track the animal took while it was at liberty. I'm currently having trouble downloading data from two different servers which I never had a problem before I recevied an upgraded computer with McAfee 8.5. When the package goes to the various servers for SST information to download, I continue to get a message in R: (...) Into my web browser, I get prompted with a file download window which quickly downloads an empty file. It looks like the data is available, I just can't seem to download it. I've deduced that it must be a problem with the new firewall features of McAfee 8.5 and was wondering if anyone else has run into a problem like this. I've tried using the exclusions options in 8.5 to exclude any R files but have come up empty handed. Still waiting on our informatics division to try and solve the problem as well but I thought I might try asking here as well. Thanks and I appreciate any assistance! Warren N. Joyce Hi Joice, In my opinion you need think about 2 things. 1- This is a problem of McAfee 8.5 not a R problem, so do you contact a McAfee support for fix this problem ? 2- I presume you use a windows but i do know your version (XP or Vista or win 95 etc) neither your R configuration (do you use R --internet2 or not?) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Division ?
On Sat, 2009-05-02 at 10:34 -0700, bogdanno wrote: It seems division with numbers bigger than 10 000 000 doesn't work 2000/21 [1] 952381 /23 [1] 2415459 Thank you Hi bogdanno First of all look this all.equal(21*2000/21,2000) [1] TRUE So de division work correctly but if I type all.equal(952381*21,2000) [1] Mean relative difference: 5e-08 It's not means R division don't work correctly if you use the command format(2000/21,digits=22) [1] 952380.952380952 So the result: 952381 is a round number not a real result of division. This occur because R print only 7 significants digits in numbers, if you test all.equal(21*952380.952380952,2000) [1] TRUE -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] omit empty cells in crosstab?
On Fri, 2009-04-24 at 13:12 -0700, sjaffe wrote: small example: a-c(1.1, 2.1, 9.1) b-cut(a,0:10) c-data.frame(b,b) d-table(c) dim(d) ##result: c(10, 10) But only 9 of the 100 cells are non-zero. If there were 10 columns, the table have 10 dimensions each of length 10, so have 10^10 elements, too much even to fit in memory Hi Steve In your only 3 cells 0 d b.1 b(0,1] (1,2] (2,3] (3,4] (4,5] (5,6] (6,7] (7,8] (8,9] (9,10] (0,1] 0 0 0 0 0 0 0 0 0 0 (1,2] 0 1 0 0 0 0 0 0 0 0 (2,3] 0 0 1 0 0 0 0 0 0 0 (3,4] 0 0 0 0 0 0 0 0 0 0 (4,5] 0 0 0 0 0 0 0 0 0 0 (5,6] 0 0 0 0 0 0 0 0 0 0 (6,7] 0 0 0 0 0 0 0 0 0 0 (7,8] 0 0 0 0 0 0 0 0 0 0 (8,9] 0 0 0 0 0 0 0 0 0 0 (9,10] 0 0 0 0 0 0 0 0 0 1 If you desire use simple code to find only cell0 use this table(interaction(c,drop=T)) (1,2].(1,2] (2,3].(2,3] (9,10].(9,10] 1 1 1 -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Dummy (factor) based on a pair of variables
On Sat, 2009-04-18 at 08:55 +0200, Serguei Kaniovski wrote: Dear All! my data is on pairs of countries, i and j, e.g.: y,i,j 1,AUT,BEL 2,AUT,GER 3,BEL,GER I would like to create a dummy (indicator) variable for use in regression (using factor?), such that it takes the value of 1 if the country is in the pair (i.e. EITHER an i-country OR an j-country). Thank you for your help, Serguei Hi Serguei, If I understand your doubt, the solution is something like this for pair i-country is AUT or j-country is BEL output ~ I(i-country==AUT|j-country==BEL) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] newbie query: simple crosstabs
On Tue, 2009-04-07 at 16:33 -0400, Donald Braman wrote: I've been playing around with various table tools, trying to construct a fairly simple cross-tab. It shouldn't be hard, but for some reason it turning out to be (for me). If I want to see how many men and how many women agree with a agree/disagree question (coded 1,0), I can do this: attach(mydata) mytable - table(male, q1.bin) # gender and a binary response variable prop.table(mytable, 1) # row percentages q1.bin male 0 1 0 0.3988 0.6012 1 0.2879 0.7121 I can repeat that for each of the items I want gender breakdowns for (q2, q3, q4 ). But what I really want is a table that shows the percentage answering yes (coded as 1) across many, many binary response items. E.g., male q1.bin q2.bin q3.bin ... 0 0.6012 0.3421 0.9871 ... 1 0.7121 0.6223 0.0198 ... I've tried various combinations of apply cbind, but to no avail. It would be easy in SPSS crosstabs, but darnit, I want to use R! Donald, The other solutions is using command CrossTable of package gmodels. -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Sequences
On Tue, 2009-04-07 at 05:16 -0700, Melissa2k9 wrote:
Hi,
I am trying to make a sequence and am using a for loop for this. I want to
start off with an initial value ie S[0]=0 then use the loop to create other
values. This is what I have so far but I just keep getting error messages.
#To calculate the culmulative sums:
s-rep(0,207)#as this is the length of the
vector I know I will have
s-as.vector(s)
s[0]-0
for (i in 1:length(lambs))# where lambs is a vector of
length 207 consisting of temperature
values
{
s[i]-s[i-1]-mean(lambs)
}
I continually get the error message:
Error in s[i] - s[i - 1] - mean(lambs) : replacement has length zero
When I merely use s[i]-i-mean(lambs) it works so there is obviously
something wrong with the s[i-1] but i cant see what. All I want is for each
S[i] to be the previous value for S - the mean.
Hi Melissa,
I think this is a index problem ...
your code is:
s-rep(0,207)
s-as.vector(s)
s[0]-0
for (i in 1:length(lambs)){
s[i]-s[i-1]-mean(lambs)
}
But try this code:
s-rep(0,207)
s-as.vector(s)
s[1]-0 # not s[0]
for (i in 2:length(lambs)){ #not 1:length(lambs)
s[i]-s[i-1]-mean(lambs)
}
The vector in R begin in 1 not in 0...
--
Bernardo Rangel Tura, M.D,MPH,Ph.D
National Institute of Cardiology
Brazil
__
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.
Re: [R] How to save the variables in R to a file
On Mon, 2009-04-06 at 11:05 +0200, Cetinyürek Aysun wrote: Dear list members, I have f.mean matrix of size 500 X 83 in R, and I want to save this into a file. I used the command save(f.mean.Rdata,file=D:/Users/Ays/Documents/Results) it saves but when I open the file in notepad it is just some characters meaningless. Thank you in advance, Kind Regards, Aysun Hi Aysun, Try write.csv something like this: write.csv(f.mean.Rdata,D:/Users/Ays/Documents/Results/f_mean.csv) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] threshold distribution
On Sun, 2009-04-05 at 01:13 -0400, jim holtman wrote: Here is what I get from using 'fitdistr' in R to fit to a lognormal. The resulting density plot from the distribution seems to be a reason match to the data. x - scan() 1: 0.80010 0.72299 0.69893 0.99597 0.89200 0.69312 0.73613 1.13559 9: 0.85009 0.85804 0.73324 1.04826 0.84002 14: 1.76330 0.71980 0.89416 0.89450 0.98670 0.83571 0.73833 0.66549 22: 0.93641 0.80418 0.95285 0.76876 0.82588 27: 1.09394 1.00195 1.14976 0.80008 1.11947 1.09484 0.81494 0.68696 35: 0.82364 0.84390 0.71402 0.80293 1.02873 40: Read 39 items plot(density(x)) library(MASS) fitdistr(x, 'lognormal') meanlogsdlog -0.134806360.19118861 ( 0.03061468) ( 0.02164785) lines(dlnorm(x, -.1348, .1911), col='red') Hi Jim I agree with your solution but my plot result not fine. I obtain same result. fitdistr(x, 'lognormal') meanlogsdlog -0.134806360.19118861 ( 0.03061468) ( 0.02164785) In plot when I use points (blue) and curve (green) the fit o lognormal and density(data) is fine but when I use lines (red) the fit is bad (in attach I send a PDF of my output) Do you know why this happen? Thanks in advance -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil P.S. my script is x - scan() 0.80010 0.72299 0.69893 0.99597 0.89200 0.69312 0.73613 1.13559 0.85009 0.85804 0.73324 1.04826 0.84002 1.76330 0.71980 0.89416 0.89450 0.98670 0.83571 0.73833 0.66549 0.93641 0.80418 0.95285 0.76876 0.82588 1.09394 1.00195 1.14976 0.80008 1.11947 1.09484 0.81494 0.68696 0.82364 0.84390 0.71402 0.80293 1.02873 require(MASS) fitdistr(x, 'lognormal') pdf(adj.pdf) plot(density(x)) lines(dlnorm(x, -.1348, .1911), col='red') points(x,dlnorm(x, -.1348, .1911), col='blue') curve(dlnorm(x, -.1348, .1911), col='green',add=T) dev.off() adj.pdf Description: Adobe PDF document __ 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.
Re: [R] Installation of all packages under Unix
On Thu, 2009-04-02 at 14:37 +0200, nboe...@uni-potsdam.de wrote: Dear List-Members, does anyone know a simple way to automatically install all R packages on a unix system to the default library path? Not from inside R, it should rather work as a shell script - job at startup. Something like R cmd install -l pkgs ### where pkgs should mean all packages (is there an option to control overwriting?) Thanks a lot, N. Boehme All packages in CRAN install.packages(new.packages(),dep=T,clean=T) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] A query about na.omit
On Wed, 2009-04-01 at 16:49 +0100, Jose Iparraguirre D'Elia wrote: Dear all, Say I have the following dataset: DF x y z [1] 1 1 1 [2] 2 2 2 [3] 3 3NA [4] 4 NA 4 [5] NA 5 5 And I want to omit all the rows which have NA, but only in columns X and Y, so that I get: x y z 1 1 1 2 2 2 3 3 NA If I use na.omit(DF), I would delete the row for which z=NA, obtaining thus x y z 1 1 1 2 2 2 But this is not what I want, of course. If I use na.omit(DF[,1:2]), then I obtain x y 1 1 2 2 3 3 which is OK for x and y columns, but I wouldn't get the corresponding values for z (ie 1 2 NA) Any suggestions about how to obtain the desired results efficiently (the actual dataset has millions of records and almost 50 columns, and I would apply the procedure on 12 of these columns)? Sincerely, Jose Luis Jose Luis Iparraguirre Senior Research Economist Economic Research Institute of Northern Ireland Hi Jose Luis, I think this script is sufficient for your problem: tab-matrix(c(1,1,1,2,2,2,3,3,NA,4,NA,4,NA,5,5),ncol=3,byrow=T) tab[!is.na(tab[,1])!is.na(tab[,2]),] -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] how to increase the limit for max.print in R
On Tue, 2009-03-31 at 14:47 +0530, pooja arora wrote: Hi All, I am using DNAcopy package in R for copy number analysis of 500K chip. The final output which I get from DNA copy is too big to be printed in a file. So I am getting an error as reached getOption(max.print) -- omitted 475569 rows Can somebody please provide me the pointers with how to increase the limit for max.print . Thanks, Pooja Hi Pooja, You must use options command, something like this options(max.print=5.5E5) For more information type? ?options -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] how to increase the limit for max.print in R
On Tue, 2009-03-31 at 15:51 +0530, pooja arora wrote: Thanks, it Worked. Do you have any idea how much is the max limit for max.print Thanks and Regards, Pooja Arora Hi Pooja, In this case max is; options(max.print=.Machine$double.xmax) In my case I have a compiled R 2.8.1 in AMD phenom using Ubuntu 8.10 AMD 64 .Machine$double.xmax [1] 1.797693e+308 -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] very fast OLS regression?
On Wed, 2009-03-25 at 22:11 +0100, Dimitris Rizopoulos wrote:
check the following options:
ols1 - function (y, x) {
coef(lm(y ~ x - 1))
}
ols2 - function (y, x) {
xy - t(x)%*%y
xxi - solve(t(x)%*%x)
b - as.vector(xxi%*%xy)
b
}
ols3 - function (y, x) {
XtX - crossprod(x)
Xty - crossprod(x, y)
solve(XtX, Xty)
}
ols4 - function (y, x) {
lm.fit(x, y)$coefficients
}
# check timings
MC - 500
N - 1
set.seed(0)
x - matrix(rnorm(N*MC), N, MC)
y - matrix(rnorm(N*MC), N, MC)
invisible({gc(); gc(); gc()})
system.time(for (mc in 1:MC) ols1(y[, mc], x[, mc]))
invisible({gc(); gc(); gc()})
system.time(for (mc in 1:MC) ols2(y[, mc], x[, mc]))
invisible({gc(); gc(); gc()})
system.time(for (mc in 1:MC) ols3(y[, mc], x[, mc]))
invisible({gc(); gc(); gc()})
system.time(for (mc in 1:MC) ols4(y[, mc], x[, mc, drop = FALSE]))
Hi Dimitris and Ivo
I think this is not a fair comparison, look this
x[8,100]-NA
system.time(for (mc in 1:MC) ols1(y[, mc], x[, mc]))
user system elapsed
8.765 0.000 8.762
system.time(for (mc in 1:MC) ols2(y[, mc], x[, mc]))
Error in solve.default(t(x) %*% x) :
system is computationally singular: reciprocal condition number = 0
Timing stopped at: 0 0 0.002
system.time(for (mc in 1:MC) ols3(y[, mc], x[, mc]))
Error in solve.default(XtX, Xty) :
system is computationally singular: reciprocal condition number = 0
Timing stopped at: 0 0 0.001
system.time(for (mc in 1:MC) ols4(y[, mc], x[, mc, drop = FALSE]))
Error in lm.fit(x, y) : NA/NaN/Inf in foreign function call (arg 1)
Timing stopped at: 0 0 0.001
So routines ols2, ols3 and ols4 only functional in fully matrix if have
one NA this functions don't run.
--
Bernardo Rangel Tura, M.D,MPH,Ph.D
National Institute of Cardiology
Brazil
__
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.
Re: [R] Search Errors
On Tue, 2009-03-24 at 11:45 -0700, CE.KA wrote: Hi R users, I have a big program So in Rgui I can t see all the execution of it Is there a way to ask R if there is Errors in my program Sincerely yours Hi Normally i use 3 functions in debug R routines: trace, browser and debug. In special cases i use a debug package ( http://cran.r-project.org/doc/Rnews/Rnews_2003-3.pdf ) Article: Debugging without (too many) tears -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] fisher.test - FEXACT error 7
On Fri, 2009-03-20 at 18:29 +, Helena Mouriño wrote: Dear all, Im having an awkward problem in R. When I write the command Fisher.test(school.data,workspace=2e+07), where school.data is the matrix corresponding to the data set, I get the error message: FEXACT error 7. LDSTP is too small for this problem. Try increasing the size of the workspace. Increasing the workspace: Fisher.test(school.data,workspace=1e+10), I get a different message, but it still doesnt work: NAs in foreign function call (arg 10) In addition: Warning message: In fisher.test(dados, workspace = 1e+10, alternative = two.sided) : NAs introduced by coercion Hi Helena, In this case you can try 3 solutions: 1- chisq.test(school.data), but pay attention if expected value of any cell is 5 2- Fisher.test(school.data,workspace=2e+07,hybrid=TRUE) from Help For larger than 2 by 2 tables and 'hybrid = TRUE', asymptotic chi-squared probabilities are only used if the 'Cochran conditions' are satisfied, that is if no cell has count zero, and more than 80% of the cells have counts at least 5. 3- Use large tables approach from Sir David Cox: Law, G. R. and Cox, D. R. and Machonochie, N. E. S. and Simpson, J. and Roman, E. and Carpenter, L. M. (2001) Large tables. Biostatistics 2(2):pp. 163-171. -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] any package for connecting berkeley db in R?
On Fri, 2009-03-20 at 17:23 -0400, Zheng, Xin (NIH) [C] wrote: Hi there, Is there any such package? I searched but found none. Thanks in advance. Xin Zheng Hi Xin Do you try package DBI? -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] numeric equality
On Wed, 2009-03-18 at 11:58 -0400, Yu, Changhong wrote: Dear all, I am totally confused by the following R output, but don't have a clue for it. a - 1 - 0.2 a == 0.8 [1] TRUE a - 1 - 0.8 a == 0.2 [1] FALSE Hi Yu, First of all read FAQ 7.31 (Why doesn't R think these numbers are equal?) Second, in this case, use all.equal a - 1 - 0.8 a == 0.2 [1] FALSE all.equal(a,0.2) [1] TRUE -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Regressão linear
On Thu, 2009-03-05 at 02:20 +, Ben Bolker wrote: Sueli Rodrigues srodrigu at esalq.usp.br writes: Olá. Tenho um arquivo que a cada 6 linhas corresponde uma amostra da qual preciso dos coeficientes da regressão linear. Como faço para que o programa distinga a cada 6 linhas como uma amostra e não calcule como um todo? Estou usando a função: model=lm(y ~ x) You're more likely to get a response if you post to the list in English (even fractured English). Based on what Google translator thinks you said (you want to perform linear regressions on 6-line subsets of a data set?), here's a starting point (assuming your data are in a data frame mydata, and have column names x and y): splitdat - split(mydata,rep(1:6,length.out=nrow(mydata)) linfits - lapply(splitdata,lm,formula=y~x) coefs - sapply(linfits,coef) or something like that. Ben Bolker Hi Ben Bolker First of all I would like to thank the kindness with my countrywoman. Second in her problem each 6 rows is a subset for a linear regression so the command is splitdat - split(mydata,rep(1:(nrow(mydata)/6),each=6)) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Does The Randvar package contain a virus(Malware) ?
On Mon, 2009-03-02 at 08:41 +0100, Matthias Kohl wrote: there were only minor changes in the latest version of RandVar (from 0.6.6 to 0.6.7). Might this be a mirror problem? Best Matthias Timthy Chang wrote: Today,I update the packages in R. but AntiVir Guard dectects the Randvar package as affected file. What happen ?? Thank you for your answer. http://www.nabble.com/file/p22281513/qq.gif Or this find a False Positive Second http://www.avira.com/en/threats/section/fulldetails/id_vir/4142/heur_html.malware.html HEUR/HTML.Malware Not a virus but a code this anti virus suspected may be a virus -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] rounding problem
On Sun, 2009-03-01 at 17:51 -0800, tedzzx wrote:
Yes, round(1.5)=2. but round(2.5)=2. I want round(2.5)=3 just like the what
the excel do. Can we change the setting or do some trick so that the
computer will work like what we usually do with respect to rounding.
My system is R 2.8.1, winXP, Intel core 2 dual . Thanks.
On Sun, 2009-03-01 at 17:51 -0800, tedzzx wrote:
Yes, round(1.5)=2. but round(2.5)=2. I want round(2.5)=3 just like
the what
the excel do. Can we change the setting or do some trick so that the
computer will work like what we usually do with respect to rounding.
If I understand your problem you need using a bias round like excell.
Do you understand this is a wrong way and this is a bug of excell, but
you need this.
Well to help you i create a function for this
round.excell-function(x){
aux-round(x,1)
aux-ifelse(
all.equal(aux-trunc(aux),.5)==TRUE,
ceiling(aux),
round(aux,0))
return(aux)
}
round.excell(1.5)
round.excell(2.5)
round.excell(2.1)
round.excell(1.1)
round.excell(2.6)
round.excell(1.6)
But I need talk to you:
1- The use of routine put a bias in your analysis
2- Is not a standard (check IEC 60559 about this)
3- Do you migrating a excell bug to R
--
Bernardo Rangel Tura, M.D,MPH,Ph.D
National Institute of Cardiology
Brazil
__
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.
Re: [R] density 1?
On Mon, 2009-03-02 at 13:27 +0100, Johannes Elias wrote: Dear R-Gurus, I wonder why 'density' values as shown in hist or plot(density(x)) are sometimes over 1. How can that be? Example hist(rnorm(1000,sd=.5),freq=FALSE) The resulting plot shows density values below 1 on the y-axis. However, hist(rnorm(1000,sd=.1),freq=FALSE) shows density values over 1. How to interpret density values over 1? Greetings, Johannes Johannes, Well density is not like probability In histogram with density the area is equal de probability in you example set.seed(123) hist(rnorm(1000,sd=.1),freq=FALSE) The interval of -0.05 and 0 have density=4 but a probability of number in this interval is 4*.05=.2 the fact set.seed(123) hist(rnorm(1000,sd=.1),freq=FALSE)$density [1] 0.0998 0.2800 0.9400 1.9800 2.6000 4. [7] 4.0400 2.9200 1.6600 0.9200 0.4400 0.1000 [13] 0.0200 set.seed(123) sum(hist(rnorm(1000,sd=.1),freq=FALSE)$density) [1] 1 So the sum of probability is 1 but the sum of density 20 -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] 1.095e+09 for integers
On Sun, 2009-02-22 at 23:42 -0500, Alexy Khrabrov wrote: I've had a very long file written out by R with write.table, with fields of time values, converted from POSIXlt as.numeric. Among 2.5 million values, very few had 6 trailing zeroes, and those were output in scientific notation as in the subject. Is this the default behavior for long integers, and how can it be turned off (with all digits for any integer field in write.table)? This is important to interoperate with other languages through such text dumps, as some do not expect scientific notation for integers, only for floats. Cheers, Alexy Alexy If I understood your problem you have something like this as.numeric(12345678912345678) [1] 1.234568e+16 In this case your solution is using format, like this format(as.numeric(12345678912345678),scientific=FALSE) [1] 12345678912345678 -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
[R] System of logistics Equations
Hi R-masters I need yours help about a problema in one of may ongoing researchers. In my research the subjects (20 in total) answer 60 questions (20 type G, 20 type S and 20 type P). Which questions is classified about 3 factor (2 level each) and the subject score with 2 scale (not integer value is possible but rare): Val range -7 to 7 and other Car range 1 to 7. This a code to fake database of research Subj-rep(1:20,each=60) Sti-rep(c(G,S,P),40) SP-rep(c(S,P),60) AG-rep(c(A,P),60) Mer-rep(c(M,NM),60) Car-round(runif(120,1,7),0) Val-round(runif(120,-7,7),0) base-data.frame(Subj,Sti,SP,AG,Mer,Car,Val) In my hypothesis: logit(Sti==G) ~ SP+AG+Mer+Car+Val+SP*Car+Mer*Val+AG*Val + errorG logit(Sti==S) ~ SP+AG+Mer+Car+Val+SP*Car+Mer*Val+AG*Val + errorS logit(Sti==P) ~ SP+AG+Mer+Car+Val+SP*Car+Mer*Val+AG*Val + errorP I test and the 3 terms of error (errorG,errorS,errorP) is correlated. So I think useful adjust a system of logistic equations to tread the 3 equations and in same time to obtain estimatives of effects and uncorrelated error terms. The systemfit package fit linear system and non-linear system but is possible adjust a logistic system in R? Thanks in advance -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
[R] System of logistics Equations-correction
Hi R-masters Sorry but i make a error in script the correct code to fake database of research Subj-rep(1:20,each=60) Sti-rep(c(G,S,P),40) SP-rep(c(S,P),each=60) AG-rep(c(A,P),60) Mer-rep(c(M,NM,M,NM),each=30) Car-round(runif(120,1,7),0) Val-round(runif(120,-7,7),0) base-data.frame(Subj,Sti,SP,AG,Mer,Car,Val) In my hypothesis: logit(Sti==G) ~ SP+AG+Mer+Car+Val+SP*Car+Mer*Val+AG*Val + errorG logit(Sti==S) ~ SP+AG+Mer+Car+Val+SP*Car+Mer*Val+AG*Val + errorS logit(Sti==P) ~ SP+AG+Mer+Car+Val+SP*Car+Mer*Val+AG*Val + errorP I test and the 3 terms of error (errorG,errorS,errorP) is correlated. So I think useful adjust a system of logistic equations to tread the 3 equations and in same time to obtain estimatives of effects and uncorrelated error terms. The systemfit package fit linear system and non-linear system but is possible adjust a logistic system in R? Thanks in advance -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Generating Numbers With Certain Distribution in R
On Wed, 2009-02-11 at 04:07 +, Ben Bolker wrote: Gundala Viswanath gundalav at gmail.com writes: Dear all, Is there a way to generate K numbers of integer (K = 10^6). The maximum value of the integer is 200,000 and minimum is 1. And the occurrences of this integer follows a lognormal distribution. - Gundala Viswanath Jakarta - Indonesia Technically speaking, I don't think this is possible since lognormal variates range from 0 to infinity and are continuous. However, perhaps something like x - rlnorm(1e6,mulog=1,sdlog=1) ## pick any parameters you like y - round((x-min(x))/diff(range(x)))*19+1 That probably doesn't do exactly what you want, but you will probably need to specify your problem more precisely (and perhaps give some context) if you want a better answer ... Ben Bolker Ben, I think your routine need a little fix x - rlnorm(1e6,meanlog=1,sdlog=1) ## pick any parameters you like y - round((x-min(x)/diff(range(x)))*19+1) What you think? -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Choosing a random number between x and y
On Mon, 2009-02-09 at 06:40 -0800, Vie wrote: Hi, Ive been trying to find a function that will allow me to pull out a number between a minimum and maximum threshold. I want a random decimal number between, for example, 0 and 0.5 or 0 and 0.7. I've been searching everywhere for a function that will allow me to do this in R, but I have yet to be successful. Any help would be much appreciated. Thanks in advance Hi Vie I don't understand if you need a only random generation or mixture random generation, so i will make the 3 examples Using runif 1- Random 10 number Retween 0 and 0.5 runif(10,0,0.5) 2 -Random 20 number Retween 0 and 0.7 runif(20,0,0.7) 3- Random 40 number of mixture two random uniforme random 1 and 2 with p(random1)= 0.3 and p(random=2) = 0.7 ifelse(runif(40).3,runif(40,0,0.7),runif(40,0,0.5)) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] How to analyse and model 2 time series, when one series needs to be differenced?
On Mon, 2009-01-26 at 08:52 +, Andreas Klein wrote: Hello. How can I analyse the cross-correlation between two time series with ccf, if one of the time series need to be differenced, so it is stationary? The two time series differ when in length and maybe ccf produces not the correct cross-correlation?! Another problem: How can I model the two time series as an VARI-process with the dse package? - So how can I handle it, that one series has to be differenced and the other series not? I hope you can give me some hints. Regards, Andreas. Hi Andreas, If I understand your problem this script solve tour question t-1:15 x-rnorm(10) y-.2-.3*t+rnorm(15) y.dif-diff(y,1) ccf(x,y.dif) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] R for Computational Neuroscience?
On Fri, 2009-01-23 at 08:53 -0400, Mike Lawrence wrote: Hi all, I've noticed that many computational neuroscience research groups use MATLAB. While it's possible that MATLAB may have some features unavailable in R, I suspect that this may instead simply be a case of costly tradition, where researchers were taught MATLAB as students and pay for it as researchers because it's all they know. I'd like to attempt to break the cycle by offering colleagues resources on using R for computational neuroscience, but I haven't been able to find anything (searched the task view, r-seek, google). Can anyone direct me to resources on using R for computational neuroscience? Input on my possibly naive assumption that R is a sufficient tool for this field would also be appreciated. Cheers, Mike Mike, I think neuroscience is a term using for a wide group of researchers. The common analysis (hypothesis test, ANOVA, regression models, etc) is perfectly made in R. But the interpretation of mri is need a packages: 1- AnalyzeFMRI -Functions for I/O, visualisation and analysis of functional Magnetic Resonance Imaging (fMRI) datasets stored in the ANALYZE or NIFTI format. 2- fmri - contains R-functions to perform an fmri analysis as described in Tabelow, K., Polzehl, J., Voss, H.U., and Spokoiny, V. Analysing fMRI experiments with structure adaptive smoothing procedures, NeuroImage, 33:55-62 (2006) 3- dti - Diffusion Weighted Imaging is a Magnetic Resonance Imaging modality, that measures diffusion of water in tissues like the human brain. The package contains R-functions to process diffusion-weighted data in the context of the diffusion tensor model (DTI). This includes the calculation of anisotropy measures and, most important, the implementation of our structural adaptive smoothing algorithm as described in K. Tabelow, J. Polzehl, V. Spokoiny, and H.U. Voss, Diffusion Tensor Imaging: Structural Adaptive Smoothing, Neuroimage 39(4), 1763-1773 (2008). -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Is it possible for R to import a SigmaPlot file?
On Thu, 2009-01-22 at 14:58 -0800, Jason Rupert wrote: I recently received a Sigmaplot file (*.jnb)from a customer and would like to know if I can input it to a data frame and then manipulate the data in R. I did a search on Google and on RSeek (www.rseek.org), but did not get any good hits. Thank for any feedback and insight you can provide. P.S. Love the flexibility of R and would love to keep using it. Just wanting to know if this is possible. Thanks again. Hi Jason, I don't know sigmaplot, but I thing is possible export sigmaplot database for other type of files. If you export jnb file to csv file is possible read a database in R. I will talk with a person work with me and tonight I send other mail with more details. -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] pchisq error
On Mon, 2009-01-19 at 09:54 +0100, Jeremy Silver wrote: Dear R experts, (...) pchisq(5.464342,1,lower.tail = FALSE) [1] 0.01940836 reproduceError(5.464342) stat = 5.464342, p = 0.019408 pchisq(5.464342,1,lower.tail = FALSE) [1] NaN Warning messages: 1: In pchisq(5.464342, 1, lower.tail = FALSE) : ** NON-convergence in pgamma()'s pd_lower_cf() f= nan. 2: In pchisq(q, df, lower.tail, log.p) : NaNs produced (...) sessionInfo() R version 2.8.0 (2008-10-20) i686-pc-linux-gnu locale: LC_CTYPE=en_US.UTF-8;LC_NUMERIC=C;LC_TIME=en_US.UTF-8;LC_COLLATE=en_US.UTF-8;LC_MONETARY=C;LC_MESSAGES=en_US.UTF-8;LC_PAPER=en_US.UTF-8;LC_NAME=C;LC_ADDRESS=C;LC_TELEPHONE=C;LC_MEASUREMENT=en_US.UTF-8;LC_IDENTIFICATION=C attached base packages: [1] stats graphics grDevices utils datasets methods base Hi Jeremy, In my computer your error is not occur. Look This: pchisq(5.464342,1,lower.tail = FALSE) [1] 0.01940836 pchisq(5.464342,1,lower.tail = FALSE) [1] 0.01940836 sessionInfo() R version 2.8.1 Patched (2009-01-16 r47630) x86_64-unknown-linux-gnu locale: LC_CTYPE=en_US.UTF-8;LC_NUMERIC=C;LC_TIME=en_US.UTF-8;LC_COLLATE=en_US.UTF-8;LC_MONETARY=C;LC_MESSAGES=en_US.UTF-8;LC_PAPER=en_US.UTF-8;LC_NAME=C;LC_ADDRESS=C;LC_TELEPHONE=C;LC_MEASUREMENT=en_US.UTF-8;LC_IDENTIFICATION=C attached base packages: [1] stats graphics grDevices utils datasets methods base Well, do you already try update your R? -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] easiest way to integrate own functions on startup
On Tue, 2009-01-20 at 01:13 +0100, Jörg Groß wrote:
Hi,
I am currently writing some own functions that I frequently need.
So, it would be perfect if I could load these functions at the
beginning of each R-session with a small command.
I tried to generate a R-package and install it that way.
But it seems that it is not so easy to add new functions to an
existing R-package.
So I am not so flexible by that.
Is there a way to just load an .R-file which contains the function-
definitions with a small command?
So that the functions are useable in the current session?
I tried load() but I get an error message...
Hi Jörg,
I do you not use .Fisrt.
Something like this:
.First - function(){
source(MyFunction1.txt)
source(MyFunction2.txt)
source(MyFunction3.txt)
...
source(MyFunctionn.txt)
}
--
Bernardo Rangel Tura, M.D,MPH,Ph.D
National Institute of Cardiology
Brazil
__
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.
Re: [R] How-To capture and handle errors with R
On Fri, 2008-12-12 at 12:51 +0100, mau...@alice.it wrote: The following error occurs every now and then by calling a function of wmTSA package: Error in `row.names-.data.frame`(`*tmp*`, value = c(1, 0)) : invalid 'row.names' length I would greatly appreciate some guidelines about how to catch such an error upon its occurrence and have it handled by my own routine rather than letting R stop the currently run script. I had a look a the R on-line documentation about errors handler and also ran some provided examples that confused my mind. The error is in textual form. I do not know whether I have to pass the whole message to the handler. When such an exception occurs the execution control is not transferred to my routine. I would like to see a working example. Thank you so much, Maura Maura you can use try to catch error for more information digit: ?try -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Logical inconsistency
Hi Kenn, Well I think your use or round isn't optimal solution. If you using round(x,1)-round(x,1) you create 2 problems First: error propagation because you make 2 round. Second: you don't using guard digits approach. The optimal use of round is using in last calculation: Look this round(8.8,1)-round(7.8,1)1 [1] TRUE round(8.8-7.8,1)1 [1] FALSE round(8.8-7.8,1)==1 [1] TRUE Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil -- Original Message --- From: Kenn Konstabel [EMAIL PROTECTED] To: emma jane [EMAIL PROTECTED] Cc: R help [EMAIL PROTECTED] Sent: Thu, 11 Dec 2008 11:53:01 +0200 Subject: Re: [R] Logical inconsistency Rounding can do no good because round(8.8,1)-round(7.8,1)1 # still TRUE round(8.8)-round(7.7)1 # FALSE What you might do is compute a-b-1 and compare it to a very small number: (8.8-7.8-1) 1e-10 # TRUE K On Wed, Dec 10, 2008 at 11:47 AM, emma jane [EMAIL PROTECTED] wrote: Thanks Greg, that does make sense. And I've solved the problem by rounding the variables before taking the difference between them. Thanks to all who replied. Emma Jane From: Greg Snow [EMAIL PROTECTED] .com.br; Wacek Kusnierczyk [EMAIL PROTECTED]; Chuck Cleland [EMAIL PROTECTED] Cc: R help [EMAIL PROTECTED] Sent: Tuesday, 9 December, 2008 16:30:08 Subject: RE: [R] Logical inconsistency Some (possibly all) of those numbers cannot be represented exactly, so there is a chance of round off error whenever you do some arithmetic, sometimes the errors cancel out, sometimes they don't. Consider: print(8.3-7.3, digits=20) [1] 1.001 print(11.3-10.3, digits=20) [1] 1 So in the first case the rounding error gives a value that is slightly greater than 1, so the greater than test returns true (if you round the result before comparing to 1, then it will return false). In the second case the uncertainties cancelled out so that you get exactly 1 which is not greater than 1 an so the comparison returns false. Hope this helps, -- Gregory (Greg) L. Snow Ph.D. Statistical Data Center Intermountain Healthcare [EMAIL PROTECTED] 801.408.8111 -Original Message- From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED] project.org] On Behalf Of emma jane Sent: Tuesday, December 09, 2008 7:02 AM To: Bernardo Rangel Tura; Wacek Kusnierczyk; Chuck Cleland Cc: R help Subject: Re: [R] Logical inconsistency Many thanks for your help, perhaps I should have set my query in context ! I'm simply calculating an indicator variable [0,1] based on the whether the difference between two measured variables is 1 or =1. I understand the FAQ about floating point arithmetic, but am still puzzled that it only apparently applies to certain elements, as follows: 8.8 - 7.8 1 TRUE 8.3 - 7.3 1 TRUE However, 10.2 - 9.2 1 FALSE 11.3 - 10.31  FALSE Emma Jane From: Bernardo Rangel Tura [EMAIL PROTECTED] To: Wacek Kusnierczyk [EMAIL PROTECTED] Cc: R help [EMAIL PROTECTED] Sent: Saturday, 6 December, 2008 10:00:48 Subject: Re: [R] Logical inconsistency On Fri, 2008-12-05 at 14:18 +0100, Wacek Kusnierczyk wrote: Berwin A Turlach wrote: Dear Emma, On Fri, 5 Dec 2008 04:23:53 -0800 (PST) Please could someone kindly explain the following inconsistencies I've discovered__when performing logical calculations in R: 8.8 - 7.8 1 TRUE 8.3 - 7.3 1 TRUE Gladly: FAQ 7.31 http://cran.at.r-project.org/doc/FAQ/R-FAQ.html#Why-doesn_0027t-R- th ink-these-numbers-are-equal_003f well, this answer the question only partially. this explains why a system with finite precision arithmetic, such as r, will fail to be logically correct in certain cases. it does not explain why r, a language said to isolate a user from the underlying implementational choices, would have to fail this way. there is, in principle, no problem in having a high-level language perform the computation in a logically consistent way. for example, bc is an arbitrary precision calculator language, and has no problem with examples as the above: bc 8.8 - 7.8 1 # 0, meaning 'no' bc 8.3 - 7.3 1 # 0, meaning 'no' bc 8.8 - 7.8 == 1 # 1, meaning 'yes' the fact that r (and many others, including matlab and sage, perhaps not mathematica) does not perform logically here is a consequence of its implementation of floating point arithmetic. the faq you were pointed to, and its referring to the goldberg's article, show that r does not successfully isolate a user from details of the lower
Re: [R] Logical inconsistency
On Thu, 2008-12-11 at 17:23 +0200, Kenn Konstabel wrote: Hi, I agree -- and my examples using round were meant as bad and dangerous examples. Using round at the last step is better and may solve the problem, but in your example ... round(8.8-7.8,1)==1 [1] TRUE ... you have to know in advance how many decimal places can possibly make a difference (is it just one? maybe it is 2? 3? 15?). round(8.8-7.8,14)==1 [1] TRUE round(8.8-7.8,15)==1 [1] FALSE ... or, equivalently, 8.8-7.8-1 1e-15 [1] TRUE 8.8-7.8-1 1e-16 [1] FALSE Best regards, Kenn Hi Kenn, If you know guard digits approach you know the answer for your question. You need the minimum numbers of significant digits and put all numbers in same number of significant digits for comparison. In all example put in this thread we need only 1 decimal place so you round for 1 decimal place, if for your job you need 3 decimal places precision you round for 3 decimal place. I think don't make sense you using 10 decimal place precision if your problem need 2 decimal places precision ... -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Logical inconsistency
On Tue, 2008-12-09 at 06:02 -0800, emma jane wrote: Many thanks for your help, perhaps I should have set my query in context ! I'm simply calculating an indicator variable [0,1] based on the whether the difference between two measured variables is 1 or =1. I understand the FAQ about floating point arithmetic, but am still puzzled that it only apparently applies to certain elements, as follows: 8.8 - 7.8 1 TRUE 8.3 - 7.3 1 TRUE However, 10.2 - 9.2 1 FALSE 11.3 - 10.31 FALSE Emma Jane Emma, This solve two forms: 1- use all.equal all.equal((10.2 - 9.2),1) [1] TRUE 2- use round round(10.2 - 9.2,0)1 [1] FALSE round(10.2 - 9.2,0)=1 [1] TRUE -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Logical inconsistency
On Fri, 2008-12-05 at 14:18 +0100, Wacek Kusnierczyk wrote: Berwin A Turlach wrote: Dear Emma, On Fri, 5 Dec 2008 04:23:53 -0800 (PST) emma jane [EMAIL PROTECTED] wrote: Please could someone kindly explain the following inconsistencies I've discovered__when performing logical calculations in R: 8.8 - 7.8 1 TRUE 8.3 - 7.3 1 TRUE Gladly: FAQ 7.31 http://cran.at.r-project.org/doc/FAQ/R-FAQ.html#Why-doesn_0027t-R-think-these-numbers-are-equal_003f well, this answer the question only partially. this explains why a system with finite precision arithmetic, such as r, will fail to be logically correct in certain cases. it does not explain why r, a language said to isolate a user from the underlying implementational choices, would have to fail this way. there is, in principle, no problem in having a high-level language perform the computation in a logically consistent way. for example, bc is an arbitrary precision calculator language, and has no problem with examples as the above: bc 8.8 - 7.8 1 # 0, meaning 'no' bc 8.3 - 7.3 1 # 0, meaning 'no' bc 8.8 - 7.8 == 1 # 1, meaning 'yes' the fact that r (and many others, including matlab and sage, perhaps not mathematica) does not perform logically here is a consequence of its implementation of floating point arithmetic. the faq you were pointed to, and its referring to the goldberg's article, show that r does not successfully isolate a user from details of the lower-level implementation. vQ Well, first of all for 8.-7.3 is not equal to 1 [for computers] 8.3-7.3-1 [1] 8.881784e-16 But if you use only one digit precision round(8.3-7.3,1)-1 [1] 0 round(8.3-7.3,1)-10 [1] FALSE round(8.3-7.3,1)==1 [1] TRUE So the problem is the code write and no the software -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] how to ignore errors
On Mon, 2008-12-01 at 16:17 +0800, Jiang Peng wrote:
Dear listers,
Hi .
I googled before I ask this question in case to avoid violating the
list law. The question was proposed by somebody else but nobody
gives a proper solution.
How to ignore errors in a R loop? I have a R source file ,namely,
test.r.
I ran a loop in a command line
for(i in 1 : 100 ) {source(test.r,echo =T ) }
and I found that when test.r have some error the loop will stop, but
test.r is just stochastically wrong,i.e., sometimes it works well if
no bad values return.
I want to know how to ignore such errors and let the loop continue.
thanks in advance.
___
Jiang Peng,
The default solution is using try command, something like this
for(i in 1 : 100 ) {
try(source(test.r,echo =T ))
}
--
Bernardo Rangel Tura, M.D,MPH,Ph.D
National Institute of Cardiology
Brazil
__
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.
Re: [R] help for code in jump diffusion
On Fri, 2008-11-21 at 22:37 -0800, subbudas wrote:
hello everyone ,
i have written some code in R for jump diffusion model.
the code generates answer as
NaN
There were 50 or more warnings (use warnings() to see the first 50)
my code is
mu-0.2
sig-0.2
S0-100
j-0.2
dt-1/252
int-0.1
i-0
while(i=1)
{
is.nan
k-rnorm(1,0,1)
theta-ifelse((k(int*dt)),1,0)
m-rnorm(1)
gam-qnorm(m,0,1)
S0-abs(S0*((1+mu*dt+sig*sqrt(dt)+ gam)- j*theta))
if(!is.nan (S0 = 0))
warning(S0 must be positive)
cat(NaN,\n)
cat(S0,\n)
i-i+(1/252)
}
the problem i am facing is i am not able to find out the reason for this NaN
output.
please help
thanks in advance.
Hi
I see two problems in your script:
1- m-rnorm(1) will produce a random number with distribution normal
mean =0 and sd=1 so m can 1 or 0
In this cases gam-qnorm(m,0,1) is NAN because m is not a probability
2- I think if(!is.nan (S0 = 0)) ... is wrong
Try this script :
mu-0.2
sig-0.2
S0-100
j-0.2
dt-1/252
int-0.1
i-0
while(i=1){
k-rnorm(1,0,1)
theta-ifelse((k(int*dt)),1,0)
# m-rnorm(1)
m-runif(1,0,1)
gam-qnorm(m,0,1)
S0-abs(S0*((1+mu*dt+sig*sqrt(dt)+ gam)- j*theta))
if(!is.nan(S0)(S0 = 0)){
warning(S0 must be positive)
cat(NaN,\n)
}
cat(S0,\n)
i-i+(1/252)
}
--
Bernardo Rangel Tura, M.D,MPH,Ph.D
National Institute of Cardiology
Brazil
__
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.
Re: [R] Mismatch in logical result?
On Fri, 2008-11-07 at 15:53 +0530, Shubha Vishwanath Karanth wrote: Hi R, I have certain checkings, which gives FALSE, but actually it is true. Why does this happen? Note that the equations that I am checking below are not even the case of recurring decimals... 1.4^2 == 1.96 [1] FALSE 1.2^3==1.728 [1] FALSE Shubha the correct answer for us is TRUE for the computer is FALSE 1.4^2-1.96 [1] -2.220446e-16 1.2^3-1.728 [1] -2.220446e-16 but if you use all.equal all.equal(1.4^2,1.96) [1] TRUE all.equal(1.2^3,1.728) [1] TRUE More details in R FAQ 7.31 -- []s Tura __ 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.
Re: [R] How to plot with different colours
On Sun, 2008-11-02 at 07:58 +, Carlos Morales wrote: Hello everyone, I'm trying to plot 3600 points and my idea is if this value is higher than 0.35 then this point must appear in green colour, if it's smaller than -0.35 then values must appear in red and if values are between -0.35 and 0.35 they must be in yellow. I'm thinking and I'm trying many things but I don't achieve it. Any idea?. Thanks so much Carlos Morales Diego Hi Carlos I think you need use a nested ifelse, something similar this x-runif(6000,-1,1) color-ifelse(x.35,green,ifelse(x -.35,red,yellow)) table(col) plot(1:6000,x,col=color) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] group sequential analysis - stopping for futility
Em Ter, 2008-10-28 às 08:30 -0700, [EMAIL PROTECTED] escreveu: Hello - I am wondering if anyone has written some R code to calculate futility bounds for a group sequential analysis of clinical trials data. The library ldbounds has a function 'bounds' which calculates the effectiveness stopping bounds for various spending functions, but it does not appear to do the same for futility (sometimes called the 'inner wedge'). Thanks for any assistance! -Alice Hi alice, I have interest in develop this code with you . Do you have a paper about this? -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] about granova library
Em Qui, 2008-10-16 às 21:11 +1030, Fernando Marmolejo Ramos escreveu: Dear all Recently the granova package was launched. I installed but after when I invoked it in R it requested for other libraries. They were downloaded and install automatically. I tried to run the example syntax of “granova.1w” and “granova.2w” but two things happened: i) either a file called “granova.rdb” wasn’t existent or ii) the GUI clashed and R shut down. Has anyone else experience this? Do the developers have an answer for this troubleshot? I’m using a Windows Vista system and I have the R version 2.7.2. Cheers, Fer Fernando I using R version 2.7.2 and Ubuntu 8.04 in my computer: granova.1w - runs fine granova.2w - don't run fine, actual only 1 of 2 graphical windows apear a plot (rgl surface) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Saving results of Kruskal Walis test
Em Qui, 2008-10-16 às 22:31 +0200, Himanshu Ardawatia escreveu: Hello, I am running Kruskal-Walis test in R. When I try to save results using write.table it gives me the following error : Error in as.data.frame.default(x[[i]], optional = TRUE, stringsAsFactors = stringsAsFactors) : cannot coerce class htest into a data.frame The overall code is as follows : data_file = read.table(~/DATA.dir/data_file.txt, header=T) attach(data_file) data_file.out - krukal.test(data_file) write.table(data_file.out, ~/DATA/results/data_file_out.txt) Error in as.data.frame.default(x[[i]], optional = TRUE, stringsAsFactors = stringsAsFactors) : cannot coerce class htest into a data.frame Results do come in data_file.out after analysis as seen below: data_file.out Kruskal-Wallis rank sum test data: value by pathway Kruskal-Wallis chi-squared = 5.6031, df = 3, p-value = 0.1326 I am wondering if I am making a mistake with using write.table (It works very well saving results from anova analysis) or is there any other way to save results in a file for future use.. Thanks Himanshu Hi Himanshu Well the output of htests is a list so data_file.out is a lista to. You don't put a list ins a data.frame so you need make this data.frame(unlist(data_file.out)) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] help about how can R compute AIC?
Em Ter, 2008-10-14 às 17:13 +0200, Arnau Mir Torres escreveu: Hello. I need to know how can R compute AIC when I study a regression model? For example, if I use these data: growth tannin 1 12 0 2 10 1 3 8 2 4 11 3 5 6 4 6 7 5 7 2 6 8 3 7 9 3 8 and I do model - lm (growth ~ tannin) AIC(model) R responses: 38.75990 I know the following formula to compute AIC: AIC= -2*log-likelihood + 2*(p+1) In my example, it would be: AIC=-2*log-likelihood + 2*2 but I don't know how R computes log-likelihood: logLik(model) 'log Lik.' -16.37995 (df=3) Arnau, LogLik= -16.37995 AIC= -2*log-likelihood + 2*(p+1) AIC=-2*-16.37995 + 2*(p+1) AIC= 32.7599+2*(p+1) # # this is very important the model have two # parameter, because sigma is a parameter to. # so # AIC= 32.7599+2*(2+1) AIC= 32.7599+6 AIC= 38.7599 -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] legend
Em Ter, 2008-10-14 às 18:44 -0700, Lavan escreveu: Hi R users, I'm trying to have the symbols for sigma[1] in my legend. the code is given below, please have a look. Thanks, lavan try legend(bottom,legend=expression(sigma[1]),...) legend(bottom,legend=paste(Sigma1=, c(0.01,0.1,0.2,0.5,1,1.5,2,4,6,9.5), sep=), fill=c(red,green,blue,black,pink,brown,purple,yellow,lightblue,orange)) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] interpreting Shapiro-Wilks test result
Em Qua, 2008-10-08 às 17:41 -0700, Halizah Basiron escreveu: Hi all, I am newbie in using R software and also doing statistical test. I want to know if my data in in normal distribution. I have 2 groups of data and I did calculate Shapiro Wilks using R software. Here is the results: Group 1: W = 0.9206, p-value = 0.01683 Group 2: W = 0.9626, p-value = 0.4694 I am not quite sure what default confidence level (CF) is used in calculating Shapiro Wilks. Else, may I choose my own CF. Let say, if I choose CF = 0.01, therefore I may approve that Group 1 and Group 2 data are normal. If not, what should I do. My next plan is to apply T Test for both groups (assuming both are normal). I really need advice from experts here. Cheers, Halizah Hi Halizah, In my opinion don't exist confidence level of test. What exist is confidence level of your experiment. If you plan your experiment with alpha probability you have 100*(1-alpha) confidence level and you p-value need minor than alpha to named significant So if you need 99% of significance the two groups is normal and you using t.test. -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Load a program at the front end
Em Qui, 2008-10-02 às 14:36 -0400, Gang Chen escreveu: I want to run a R program, prog.R, interactively. My question is, is there a way I can start prog.R on the shell terminal when invoking R, instead of using source() inside R? TIA, Gang Hi Gang I my system just only type: R --no-save prog.R platform x86_64-unknown-linux-gnu arch x86_64 os linux-gnu system x86_64, linux-gnu status Patched major 2 minor 7.2 year 2008 month 09 day11 svn rev46532 language R version.string R version 2.7.2 Patched (2008-09-11 r46532) -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] cannot open connection: Authorization Required
Em Sex, 2008-10-03 às 11:22 -0700, Spencer Graves escreveu:
Hi, All:
Is there a way in R to access a file / web site that requires
permission?
Consider for example the following:
readLines('http://www.r-project.org/', 4)
[1] !DOCTYPE HTML PUBLIC \-//W3C//DTD HTML 4.01 Transitional//EN\
[2] html
[3] head
[4] titleThe R Project for Statistical Computing/title
readLines(URL) # URL = web address, which I can see via manual access.
Error in file(con, r) : cannot open the connection
In addition: Warning message:
In file(con, r) :
cannot open: HTTP status was '401 Authorization Required'
Thanks,
Spencer Graves
Spencer, this is a local problem in my system
readLines('http://www.r-project.org/')
[1] !DOCTYPE HTML PUBLIC \-//W3C//DTD HTML 4.01 Transitional//EN\
[2] html
[3] head
[4] titleThe R Project for Statistical Computing/title
[5] link rel=\icon\ href=\favicon.ico\ type=\image/x-icon\
[6] link rel=\shortcut icon\ href=\favicon.ico\
type=\image/x-icon\
[7] link rel=\stylesheet\ type=\text/css\ href=\R.css\
[8] /head
[9]
[10] FRAMESET cols=\1*, 4*\ border=0
[11] FRAMESET rows=\120, 1*\
[12] FRAME src=\logo.html\ name=\logo\ frameborder=0
[13] FRAME src=\navbar.html\ name=\contents\ frameborder=0
[14] /FRAMESET
[15] FRAME src=\main.shtml\ name=\banner\ frameborder=0
[16] noframes
[17] h1The R Project for Statistical Computing/h1
[18]
[19] Your browser seems not to support frames,
[20] here is the A href=\navbar.html\contents page/A of the R
Project's
[21] website.
[22] /noframes
[23] /FRAMESET
[24]
[25]
[26]
I run R.2.7.2 in Ubuntu AMD 64 machine
--
Bernardo Rangel Tura, M.D,MPH,Ph.D
National Institute of Cardiology
Brazil
__
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.
Re: [R] Logistic regression problem
Em Ter, 2008-09-30 às 18:56 -0500, Frank E Harrell Jr escreveu: Bernardo Rangel Tura wrote: Em Sáb, 2008-09-27 às 10:51 -0700, milicic.marko escreveu: I have a huge data set with thousands of variable and one binary variable. I know that most of the variables are correlated and are not good predictors... but... It is very hard to start modeling with such a huge dataset. What would be your suggestion. How to make a first cut... how to eliminate most of the variables but not to ignore potential interactions... for example, maybe variable A is not good predictor and variable B is not good predictor either, but maybe A and B together are good predictor... Any suggestion is welcomed milicic.marko I think do you start with a rpart(binary variable~.) This show you a set of variables to start a model and the start set to curoff for continous variables I cannot imagine a worse way to formulate a regression model. Reasons include 1. Results of recursive partitioning are not trustworthy unless the sample size exceeds 50,000 or the signal to noise ratio is extremely high. 2. The type I error of tests from the final regression model will be extraordinarily inflated. 3. False interactions will appear in the model. 4. The cutoffs so chosen will not replicate and in effect assume that covariate effects are discontinuous and piecewise flat. The use of cutoffs results in a huge loss of information and power and makes the analysis arbitrary and impossible to interpret (e.g., a high covariate value:low covariate value odds ratio or mean difference is a complex function of all the covariate values in the sample). 5. The model will not validate in new data. Professor Frank, Thank you for your explain. Well, if my first idea is wrong what is your opinion on the following approach? 1- Make PCA with data excluding the binary variable 2- Put de principal components in logistic model 3- After revert principal componentes in variable (only if is interesting for milicic.marko) If this approach is wrong too what is your approach? -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] birthday problem (factorial limit)
Em Dom, 2008-09-28 às 19:43 +0200, Jörg Groß escreveu: Hi, I tried to calculate the formula for the birthday problem (the probability that at least two people out of a group of n people share the same birthday) But the factorial-function allows me only to calculate factorials up to 170. So is there a way to push that limit? to solve this formula: (factorial(365) / factorial((365-23))) / (365^23) (n=23) Log experession n-23 exp(sum(log(1:365))-sum(log(1:(365-n)))-n*log(365)) [1] 0.4927028 -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Logistic regression problem
Em Sáb, 2008-09-27 às 10:51 -0700, milicic.marko escreveu: I have a huge data set with thousands of variable and one binary variable. I know that most of the variables are correlated and are not good predictors... but... It is very hard to start modeling with such a huge dataset. What would be your suggestion. How to make a first cut... how to eliminate most of the variables but not to ignore potential interactions... for example, maybe variable A is not good predictor and variable B is not good predictor either, but maybe A and B together are good predictor... Any suggestion is welcomed milicic.marko I think do you start with a rpart(binary variable~.) This show you a set of variables to start a model and the start set to curoff for continous variables -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] cannot allocate memory
Em Ter, 2008-09-23 às 21:42 -0400, DumpsterBear escreveu: I am getting Error: cannot allocate vector of size 197 MB. I know that similar problems were discussed a lot already, but I didn't find any satisfactory answers so far! Details: *** I have XP (32bit) with 4GB ram. At the time when the problem appeared I had 1.5GB of available physical memory. *** I increased R memory limit to 3GB via memory.limit(3000) *** I did gs() and got used (Mb) gc trigger (Mb) max used (Mb) Ncells147534 4.0 407500 10.9407500 10.9 Vcells 104939449 800.7 186388073 1422.1 185874684 1418.2 The garbage collection didn't help. Any ideas? Many thanks in advance! Adam, First, is possible 32bit XP use all your 4Gb? Second, I think you say gc whem say gs, so in my computer (Ubuntu 64bit with 4Gb): gc() used (Mb) gc trigger (Mb) max used (Mb) Ncells 188975 10.1 407500 21.8 35 18.7 Vcells 169133 1.3 786432 6.0 786378 6.0 gc(reset=T) used (Mb) gc trigger (Mb) max used (Mb) Ncells 188951 10.1 407500 21.8 188951 10.1 Vcells 168893 1.3 786432 6.0 168893 1.3 If you read gc help: reset: logical; if 'TRUE' the values for maximum space used are reset to the current values. Other issue is options for rgui command. Have a option --max-mem-size that you modify to expand you RAM avaiable -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
Re: [R] Confidence Interval
Em Ter, 2008-08-19 às 23:25 -0300, Raphael Saldanha escreveu: Hi! With the following script, I'm trying to make a demonstration of a Confidence Interval, but I'm observing some differences on tails. Raphael, If you make demonstration of Confidence Interval why you don't use ci.examp of TeachingDemos package? It's a very good example -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil __ 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.
