Re: [R] How to use Sys.time() while writing a csv file name
Dear Mr Newmiller and Mr Oettli,
Thanks a lot for your valuable guidance. Task is done. Thanks again.
Regards
Vincy
--- On Wed, 7/4/12, Jeff Newmiller jdnew...@dcn.davis.ca.us wrote:
From: Jeff Newmiller jdnew...@dcn.davis.ca.us
Subject: Re: [R] How to use Sys.time() while writing a csv file name
To: Vincy Pyne vincy_p...@yahoo.ca, r-help@r-project.org
Received: Wednesday, July 4, 2012, 5:38 AM
You forgot to follow the posting guide and tell us what operating system you
are using (sessionInfo), but I am going to guess that you are on Windows where
the colon (:) is an illegal symbol in filenames. Try formatting the time
explicitly in the conversion to character using the format string definitions
found in ?strptime in a format that doesn't include colons.
---
Jeff Newmiller            The    . Â
   . Go Live...
DCN:jdnew...@dcn.davis.ca.us    Basics: ##.#.     ##.#. Live
Go...
                   Live:   OO#.. Dead:
OO#.. Playing
Research Engineer (Solar/Batteries      O.O#.     #.O#.Â
with
/Software/Embedded Controllers)         .OO#.     .OO#.Â
rocks...1k
---
Sent from my phone. Please excuse my brevity.
Dear R helpers,
I am using Beta distribution to generate the random no.s (recovery
rates in my example). However, each time I need to save these random
no.s in a csv format. To distinguish different csv files, one way I
thought was use of Sys.time in the file name. My code is as follows -
# My code
rr = rbeta(25, 6.14, 8.12)
lgd = 1 - mean(rr)
write.csv(data.frame(recovery_rates = rr), file =
paste(recovery_rates_at_, Sys.time(), .csv, sep = ), row.names =
FALSE)
However, I get following error -
Error in file(file, ifelse(append, a, w)) : � cannot open the
connection
In addition: Warning message: In file(file, ifelse(append, a, w)) :
cannot open file 'recovery_rates_at_2012-07-04 1:14:05.csv': Invalid
argument
If instead of Sys.time, I use some other variable e.g. lgd as
write.csv(data.frame(recovery_rates = rr), paste('rates_',lgd,'.csv',
sep = ), row.names = FALSE)
I am able to store these simulated recovery rates in different files.
But I need to use Sys.time in my csv file name. (or is there any other
way of writing these csv files so that files don't get over-written).
Kindly guide.
Regards and thanking in advance
Vincy
   [[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.
[[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.
[R] problem loading siar
Hi, I have a problem while loading the siar program in R. When I am loading siar, system does not load convexhull. On the screen I have seen such writings. The following object(s) are masked from package:spatstat: convexhull How can I load the convexhull, how can I unmask from this package? I will be appreciated if you give advice about this. Best Sukran [[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.
Re: [R] problem loading siar
Hello, It is not what happens. Function convexhull exists in both siar and spatstat packages. As you already load spatstat, when you are loading siar, the convexhull in spatstat is masked by the one in siar. Thus, when you will run convexhull function, it will be the one from the siar package. Regards Le 04/07/2012 15:24, Sukran yalcin ozdilek a écrit : Hi, I have a problem while loading the siar program in R. When I am loading siar, system does not load convexhull. On the screen I have seen such writings. The following object(s) are masked from ‘package:spatstat’: convexhull How can I load the convexhull, how can I unmask from this package? I will be appreciated if you give advice about this. Best Sukran [[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. __ 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 filled.contour() -
Dear all, I can't figure out a way to have more than one plot using filled.contour() in a single plate. I tried to use layout() or par(), but the way filled.contour() is written seems to override those commands. Any suggestions would be really appreciated. Jonathan [[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.
Re: [R] Data manipulation with aggregate
Hi,
Try this:
myData = data.frame(Name = c('a', 'a', 'b', 'b'), length = c(1,2,3,4), type=
c('x','x','y','z'))
z-aggregate(length~Name,myData,mean)
z1-aggregate(length~type,myData,mean)
merge(z,merge(z,z1),all=TRUE)
Name length type
1 a 1.5 x
2 b 3.5 NA
A.K.
- Original Message -
From: Filoche pmassico...@hotmail.com
To: r-help@r-project.org
Cc:
Sent: Tuesday, July 3, 2012 12:04 PM
Subject: [R] Data manipulation with aggregate
Hi everyone.
I have these data :
myData = data.frame(Name = c('a', 'a', 'b', 'b'), length = c(1,2,3,4), type
= c('x','x','y','z'))
which gives me:
Name length type
1 a 1 x
2 a 2 x
3 b 3 y
4 b 4 z
I would group (mean) this DF using 'Name' as grouping factor. However, I
have a field ('type') which is a string. I would like to use the unique
value of this field when possible (i.e. when all the 'type' values are the
same for each group) or replace with NA when 'type' has multiple values.
In fact, I would like to obtain this:
Name length type
1 a 1.5 x
2 b 3.5 NA
For instance, I was using this command:
aggregate(list(myData$length, myData$type), list(myData$Name), FUN = mean)
But it can't deal with string data.
I hope I have been clear enough.
With regards,
Phil
--
View this message in context:
http://r.789695.n4.nabble.com/Data-manipulation-with-aggregate-tp4635298.html
Sent from the R help mailing list archive at Nabble.com.
__
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@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] Please help
Dear All,
I am a research student in environment. I have only little
programming knowledge. I am currently doing the last project about rainfall
impact on ground water quality in an area. It happens that I have to use R to
read rainfall data (3 dimension) from ASC file (*.asc), and then write them
into one NCDF file (*.nc).
I have been working very hard on study R, but I still can not fix
the problem. Could someone please as kind as point out that what the problems
are in my R script?
Firstly, this is an example of data in asc file:
NCOLS 241
NROWS 160
XLLCORNER140.00012207031
YLLCORNER -39.
CELLSIZE0.50E-01
NODATA_VALUE -99.0
166.30 160.87 155.23 149.33 143.83 138.52 133.29
128.34 123.76 119.21
115.06 110.90 107.22 103.69 100.40 97.29 94.58
92.15 90.00 87.91
86.20 84.57 83.22 81.94 81.11 80.38 79.37
78.73 79.70 79.62
---
And then this is the script I wrote:
setwd(E:/grid)
#defining dimension
x=dim.def.ncdf(Lon,degreesE,140.0251:146.6751)
y=dim.def.ncdf(Lat,degreesN,(-31.025):(-38.975))
t=dim.def.ncdf(Time,1968-01,1:12,unlim=TRUE)
#setup variable
varmr=var.def.ncdf(mr,mm,list(x,y,t),-99.00,
longname=monthly rainfall)
#create ncdf file
ncnew=create.ncdf(rainfall.nc, varmr)
#read input
files=list.files(pattern=.asc)
mrain=matrix(0:0,0,3)
for(i in files)
{rainfall=data.frame(readAsciiGrid(i))
mrain=cbind(mrain,rainfall)
}
put.var.ncdf(ncnew, mrain)
close.ncdf(ncnew)
---
[[elided Hotmail spam]]
Many thanks,
Jun
[[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.
Re: [R] Date
Hi, Try this: Year-c(01/2000,02/2000,03/2000) #If you want to convert it directly to month/year library(zoo) as.yearmon(Year,format=%m/%Y) [1] Jan 2000 Feb 2000 Mar 2000 #As your intention is to have DD/MM/ format, Year1-paste(01/,Year,sep=) Year1 [1] 01/01/2000 01/02/2000 01/03/2000 Year-as.Date(Year1, format=%m/%d/%Y) Year [1] 2000-01-01 2000-01-02 2000-01-03 dat1-data.frame(Year1,Stock_Prices=1:3) dat1 Year1 Stock_Prices 1 01/01/2000 1 2 01/02/2000 2 3 01/03/2000 3 dat2-data.frame(Year,Stock_Prices=1:3) dat2 Year Stock_Prices 1 2000-01-01 1 2 2000-01-02 2 3 2000-01-03 3 A.K. - Original Message - From: Akhil dua akhil.dua...@gmail.com To: r-help@r-project.org Cc: Sent: Wednesday, July 4, 2012 12:13 AM Subject: [R] Date Hi I have monthly data and the dates are in MM/YY Format I need to convert them into DD/MM/YY format by pasting 01 in place of DD to all the observations in my Year Column ex: Year Stock Prices 01/2000 1 02/2000 2 03/2000 3 I need to convert them to Year Stock Prices 01/01/2000 1 01/02/2000 2 01/03/2000 3 [[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. __ 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 Sys.time() while writing a csv file name
Hello, It seems that there is a problem with :. If you only need the date, you can use as.Date(Sys.time()) instead of the complete form. If you need the time, you can try the following commands which change the : into - : newsystime--format(Sys.time(),%Y-%m-%d-%H-%M-%S) write.csv(data.frame(recovery_rates = rr), file = paste(recovery_rates_at_, newsystime, .csv, sep = ), row.names = FALSE) Have a good day, Ptit Bleu. -- View this message in context: http://r.789695.n4.nabble.com/How-to-use-Sys-time-while-writing-a-csv-file-name-tp4635358p4635363.html Sent from the R help mailing list archive at Nabble.com. __ 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] Printing from R Console in colour
Hi,I want to be able to print in colour from the R console i.e. commands (in navy) with results (in red) as on the console. From the console if I click on File - Print, the commands with results print on my printer but only in monochrome, not in colour. Similarly, if I Edit - Select All, Edit - Copy --- and paste into Word, the commands and results are only in monochrome, not in colour. How do I enable R to print commands and results in colour? Thanks [[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.
Re: [R] (no subject)
Hello,
Or maybe to avoid the typo (Market is column 5), use variables names, in
something like
myData - read.table(text=
Date Stock1 Stock2 Stock3Market
01/01/2000 1 2 3 4
01/02/2000 5 6 7 8
01/03/2000 1 2 3 4
01/04/2000 5 6 7 8
, header=TRUE, stringsAsFactors=FALSE)
myData$Date - as.Date(myData$Date, fortmat=%m/%d/%Y)
myData
# Avoid typos
stocks - grep(Stock, names(myData))
models - lapply(myData[, stocks], function(x) lm(x ~ myData$Market))
# Do whatever you want with results
lapply(models, summary)
Hope this helps,
Rui Barradas
Em 04-07-2012 06:29, Hasan Diwan escreveu:
On 3 July 2012 22:03, Akhil dua akhil.dua...@gmail.com wrote:
and I need to run a seperate regression of every stock on market
so I want to write a for loop so that I wont have to write codes again
and again to run the regression...
1. Do give a subject line -- a blank one is commonly used by a virus.
2. In R/S+/most functional languages, you do not want to write a for
loop. Use apply (and friends) instead.
my data is in the format given below
Date Stock1 Stock2 Stock3Market
01/01/2000 1 2 3 4
01/02/2000 5 6 7 8
01/03/2000 1 2 3 4
01/04/2000 5 6 7 8
For example, if you wanted to know the stocks share of the total market as
a fraction, you'd use something like:
sapply(myData[,c(2:4)], function(x) {
return(as.numeric(x)/as.numeric(myData[,4])) })
Hope that helps.. -- H
__
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] how to check convergence of arima model
I have already fitted several models using R code; arima(rates,c(p,d,q)) As I heard, best model produce the smallest AIC value, but maximum likelihood estimation procedure optimizer should converge. How to check whether maximum likelihood estimation procedure optimizer has converged or not? [[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.
Re: [R] problem loading siar
Hi Hello, It is not what happens. Function convexhull exists in both siar and spatstat packages. As you already load spatstat, when you are loading siar, the convexhull in spatstat is masked by the one in siar. Thus, when you will run convexhull function, it will be the one from the siar package. AFAIK there is an option to use both functions, you just need to specify from which package you want it. I do not use it but I believe there is a mention in docs and it was discussed before in help list too. probably spatstat:::convexhull() Regards Petr Regards Le 04/07/2012 15:24, Sukran yalcin ozdilek a écrit : Hi, I have a problem while loading the siar program in R. When I am loading siar, system does not load convexhull. On the screen I have seen such writings. The following object(s) are masked from ‘package:spatstat’: convexhull How can I load the convexhull, how can I unmask from this package? I will be appreciated if you give advice about this. Best Sukran [[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. __ 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@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] ggplot2: legend
Hi I do not have direct answer. You shall probably search ggplot2 web. Searching legend gave me about sixty results from which you probably could learn how to modify legend(s) according to your wish. e.g. http://had.co.nz/ggplot2/docs/opts.html Regards Petr Dear all, I produced the following graph with ggplot which is almost fine, yet I don't like that the legend for Means and Observations includes a line, though no line is used in the plot for those two (the line for Overall Mean on the other hand is wanted): library(ggplot2) ddf - data.frame(x = factor(rep(LETTERS[1:2], 5)), y = rnorm(10)) p - ggplot(ddf, aes(x = x, y = y)) p + geom_point(aes(colour=Observations, shape=Observations)) + stat_summary(aes(colour = Means, shape =Means), fun.y = mean, fun.ymin = min, fun.ymax = max, geom = point) + geom_hline(aes(yintercept = mean(y), linetype = Overall Mean), show_guide = TRUE) I tried to map the linetype in geom_point (and stat_summary) to NULL and I tried to set it to blank but neither worked. Is it furthermore possible to combine the two legends? My preferred plot would have two symbols with different colours and shapes for Means and Observations (but no line) and directly below that a line for Overall Mean (that is all the three items in one legend)? For that I tried to assign the same names to the legends but this did not work either. So any help would be highly appreciated. Kind Regards, Thorn Thaler Mathematician Applied Mathematics Nestec Ltd, Nestlé Research Center PO Box 44 CH-1000 Lausanne 26 Phone: +41 21 785 8220 Fax: +41 21 785 9486 __ 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@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 challenging question: merging excel files under a specific pattern
Hi Well, this is help list for R not for Excel, maybe you shall contact Microsoft guys. I believe that probably easiest would be to make a simple macro in Excel. If you want to do merging in R you shall go through help pages for read.xls, merge, cbind, rbind and R data import/export manual. Regards Petr Dear all, I have an excel file that contains 6 sheets 1,2,3,4,5,6 The analysis is repeated every 3 sheets Sheets 1, 2, 3: I want to add (horizontally) the data contained in the matrix : sheet2 (5:end,3:end ) of *Sheet2 * to the sheet3 such that the first element of the matrix *sheet2 (5:end,3:end ) * to occupy the location/cell sheet3(5,end+1 ) of sheet3. Say, that the output from this merging is sheetA. Then I want to add horizontally the data contained in the matrix : Sheet1 (5:end,3:end ) of Sheet 1 to the sheetA * such that the first element of the matrix *sheet1 (5: end,3:end ) to occupy the location/cell sheetA(5,end+1 ) of sheetA. As you can see 1)I add sheet2 (5:end,3:end ) * to *sheet3 at location *sheet3(5,end+1) * 2) then I add sheet1 (5:end,3:end ) to the output sheetA that results from the merging of sheets 2 and 3 at location sheetA((5,end+1). 3) The output is named ,say, sheetAA Similarly analysis holds for the other block of sheets 4,5,6. That is, Sheets 4, 5, 6: 1)I add sheet5 (5:end,3:end ) to sheet6 at location sheet6(5,end+1) 2) then I add sheet4 (5:end,3:end ) to the output sheetB that results from the merging of sheets 5 and 6 at location sheetB((5,end+1). 3) The output is named, say, sheetBB In my case I have a large sequence of sheets that I have to merge in this way. That is, 1,2,3,4,5,6,7,8,9,10,11,12,13,… But the logic is the same as described above. Is there any “easy” way to do that kind of merging? . Because when you have 13*3 =39 sheets is a bit tedious to do that merging manually. thanks __ 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@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 check convergence of arima model
Hello, Sorry, but do you read the answers to your posts? Inline. Em 04-07-2012 08:02, Sajeeka Nanayakkara escreveu: I have already fitted several models using R code; arima(rates,c(p,d,q)) And I have already answered to a question starting like this yesterday. In the mean time, the subject line changed. As I heard, best model produce the smallest AIC value, but maximum likelihood estimation procedure optimizer should converge. How to check whether maximum likelihood estimation procedure optimizer has converged or not? Yes, it was this question, the subject line was 'question'... ... And the answer was: read the manual, that I quoted, by the way. It now changed to: read the manual, period. Rui Barradas [[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. __ 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 filled.contour() -
Take a look at the code for filled.contour(). You'll find a line beginning .Internal(filledcontour( You can adapt this line and the lines around it to achieve what you want. Good luck Bob Kinley -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of Jonathan Hughes Sent: 04 July 2012 02:01 To: r-help@r-project.org Subject: [R] help with filled.contour() - Dear all, I can't figure out a way to have more than one plot using filled.contour() in a single plate. I tried to use layout() or par(), but the way filled.contour() is written seems to override those commands. Any suggestions would be really appreciated. Jonathan [[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. __ 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 check convergence of arima model
Hello, Inline. Em 04-07-2012 09:35, Sajeeka Nanayakkara escreveu: Hi, Sorry, since I didn't see the earlier message I resend it. I read the help page that you mentioned. But the problem is for all models, code is zero. According to that, all models were converged. Considering AIC value the best model is selected. No, arima() does not select models by AIC. That is the default behavior of ar(); arima() does NOT select models, it selects, using optim, values for the parameters of a specified model. You must choose the orders yourself. Hope this helps, Rui Barradas Is that correct procedure? The R code which was used is; model1-arima(rates,c(1,1,1)) model1 model1$code [1] 0 Sajeeka Nanayakkara *From:* Rui Barradas ruipbarra...@sapo.pt *To:* Sajeeka Nanayakkara nsaje...@yahoo.com *Cc:* r-help@r-project.org *Sent:* Wednesday, July 4, 2012 2:01 PM *Subject:* Re: [R] how to check convergence of arima model Hello, Sorry, but do you read the answers to your posts? Inline. Em 04-07-2012 08:02, Sajeeka Nanayakkara escreveu: I have already fitted several models using R code; arima(rates,c(p,d,q)) And I have already answered to a question starting like this yesterday. In the mean time, the subject line changed. As I heard, best model produce the smallest AIC value, but maximum likelihood estimation procedure optimizer should converge. How to check whether maximum likelihood estimation procedure optimizer has converged or not? Yes, it was this question, the subject line was 'question'... ... And the answer was: read the manual, that I quoted, by the way. It now changed to: read the manual, period. Rui Barradas [[alternative HTML version deleted]] __ R-help@r-project.org mailto: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@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 check convergence of arima model
Hi, I need to predict exchange rates using time series. So according to ACF and PACF knowledge, mixed model (ARIMA) be the appropriate and now, I need to find order of the model (p,d,q). So, several models were fitted to select the suitable model using arima(). Could you please tell me the procedure of selecting the correct order as I don't have enough time to search? Thank you. Sajeeka Nanayakkara From: Rui Barradas ruipbarra...@sapo.pt Cc: r-help@r-project.org Sent: Wednesday, July 4, 2012 2:58 PM Subject: Re: [R] how to check convergence of arima model Hello, Inline. Em 04-07-2012 09:35, Sajeeka Nanayakkara escreveu: Hi, Sorry, since I didn't see the earlier message I resend it. I read the help page that you mentioned. But the problem is for all models, code is zero. According to that, all models were converged. Considering AIC value the best model is selected. No, arima() does not select models by AIC. That is the default behavior of ar(); arima() does NOT select models, it selects, using optim, values for the parameters of a specified model. You must choose the orders yourself. Hope this helps, Rui Barradas Is that correct procedure? The R code which was used is; model1-arima(rates,c(1,1,1)) model1 model1$code [1] 0 Sajeeka Nanayakkara *From:* Rui Barradas ruipbarra...@sapo.pt *Cc:* r-help@r-project.org *Sent:* Wednesday, July 4, 2012 2:01 PM *Subject:* Re: [R] how to check convergence of arima model Hello, Sorry, but do you read the answers to your posts? Inline. Em 04-07-2012 08:02, Sajeeka Nanayakkara escreveu: I have already fitted several models using R code; arima(rates,c(p,d,q)) And I have already answered to a question starting like this yesterday. In the mean time, the subject line changed. As I heard, best model produce the smallest AIC value, but maximum likelihood estimation procedure optimizer should converge. How to check whether maximum likelihood estimation procedure optimizer has converged or not? Yes, it was this question, the subject line was 'question'... ... And the answer was: read the manual, that I quoted, by the way. It now changed to: read the manual, period. Rui Barradas [[alternative HTML version deleted]] __ R-help@r-project.org mailto: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. [[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.
Re: [R] how to check convergence of arima model
Hello, Put the fitted models in a list and then use lapply with AIC(). Like this set.seed(1) x - 1:100 + sqrt(1:100 + runif(100)) + rnorm(100) models - list() models[[1]] - arima(diff(x), order = c(1, 0, 0)) # Just to show models[[2]] - arima(diff(x), order = c(1, 0, 1)) # several models[[3]] - arima(diff(x), order = c(2, 0, 2)) # models models[[4]] - arima(diff(x), order = c(2, 0, 3)) lapply(models, AIC) Or run each model at a time through AIC(), whichever suits better. Hope this helps Rui Barradas Em 04-07-2012 10:22, Sajeeka Nanayakkara escreveu: Hi, I need to predict exchange rates using time series. So according to ACF and PACF knowledge, mixed model (ARIMA) be the appropriate and now, I need to find order of the model (p,d,q). So, several models were fitted to select the suitable model using arima(). Could you please tell me the procedure of selecting the correct order as I don't have enough time to search? Thank you. Sajeeka Nanayakkara *From:* Rui Barradas ruipbarra...@sapo.pt *To:* Sajeeka Nanayakkara nsaje...@yahoo.com *Cc:* r-help@r-project.org *Sent:* Wednesday, July 4, 2012 2:58 PM *Subject:* Re: [R] how to check convergence of arima model Hello, Inline. Em 04-07-2012 09:35, Sajeeka Nanayakkara escreveu: Hi, Sorry, since I didn't see the earlier message I resend it. I read the help page that you mentioned. But the problem is for all models, code is zero. According to that, all models were converged. Considering AIC value the best model is selected. No, arima() does not select models by AIC. That is the default behavior of ar(); arima() does NOT select models, it selects, using optim, values for the parameters of a specified model. You must choose the orders yourself. Hope this helps, Rui Barradas Is that correct procedure? The R code which was used is; model1-arima(rates,c(1,1,1)) model1 model1$code [1] 0 Sajeeka Nanayakkara *From:* Rui Barradas ruipbarra...@sapo.pt mailto:ruipbarra...@sapo.pt *To:* Sajeeka Nanayakkara nsaje...@yahoo.com mailto:nsaje...@yahoo.com *Cc:* r-help@r-project.org mailto:r-help@r-project.org *Sent:* Wednesday, July 4, 2012 2:01 PM *Subject:* Re: [R] how to check convergence of arima model Hello, Sorry, but do you read the answers to your posts? Inline. Em 04-07-2012 08:02, Sajeeka Nanayakkara escreveu: I have already fitted several models using R code; arima(rates,c(p,d,q)) And I have already answered to a question starting like this yesterday. In the mean time, the subject line changed. As I heard, best model produce the smallest AIC value, but maximum likelihood estimation procedure optimizer should converge. How to check whether maximum likelihood estimation procedure optimizer has converged or not? Yes, it was this question, the subject line was 'question'... ... And the answer was: read the manual, that I quoted, by the way. It now changed to: read the manual, period. Rui Barradas [[alternative HTML version deleted]] __ R-help@r-project.org mailto:R-help@r-project.org mailto:R-help@r-project.org mailto: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@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] Please help
Hello,
Following lines are wrong:
x=dim.def.ncdf(Lon,degreesE,140.0251:146.6751)
y=dim.def.ncdf(Lat,degreesN,(-31.025):(-38.975))
You have 241 longitudes and 160 latitudes.
lon - seq(from=140.0251, to=146.6751, length.out=241)
lat - seq(from=-38.975, to=-31.025, length.out=160)
Regards
- Mail original -
De : Jun Chen chensh...@hotmail.com
À : r-help@r-project.org
Cc :
Envoyé le : Mercredi 4 juillet 2012 10h52
Objet : [R] Please help
Dear All,
I am a research student in environment. I have only little
programming knowledge. I am currently doing the last project about rainfall
impact on ground water quality in an area. It happens that I have to use R to
read rainfall data (3 dimension) from ASC file (*.asc), and then write them
into one NCDF file (*.nc).
I have been working very hard on study R, but I still can not fix
the problem. Could someone please as kind as point out that what the problems
are in my R script?
Firstly, this is an example of data in asc file:
NCOLS 241
NROWS 160
XLLCORNER 140.00012207031
YLLCORNER -39.
CELLSIZE 0.50E-01
NODATA_VALUE -99.0
166.30 160.87 155.23 149.33 143.83 138.52 133.29
128.34 123.76 119.21
115.06 110.90 107.22 103.69 100.40 97.29 94.58
92.15 90.00 87.91
86.20 84.57 83.22 81.94 81.11 80.38 79.37
78.73 79.70 79.62
---
And then this is the script I wrote:
setwd(E:/grid)
#defining dimension
x=dim.def.ncdf(Lon,degreesE,140.0251:146.6751)
y=dim.def.ncdf(Lat,degreesN,(-31.025):(-38.975))
t=dim.def.ncdf(Time,1968-01,1:12,unlim=TRUE)
#setup variable
varmr=var.def.ncdf(mr,mm,list(x,y,t),-99.00,
longname=monthly rainfall)
#create ncdf file
ncnew=create.ncdf(rainfall.nc, varmr)
#read input
files=list.files(pattern=.asc)
mrain=matrix(0:0,0,3)
for(i in files)
{rainfall=data.frame(readAsciiGrid(i))
mrain=cbind(mrain,rainfall)
}
put.var.ncdf(ncnew, mrain)
close.ncdf(ncnew)
---
[[elided Hotmail spam]]
Many thanks,
Jun
[[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.
__
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] sqlSave()
Julien Moeys: ok, I see what is the problem, Your example does not work because MS Access is trying to update values in your table according to the ID you provide So when you provide for instance ID = 1, MS Access will look in the table for an existing record(s) that have an ID of 1, and replace the existing value by the new one When you provide ID = 16, MS Access will look in the table for an existing record(s) that have an ID of 16, but will not find it, thus the error So you should: - use sqlSave() for the records that have an ID that is *not yet* in the database table (in your example below, all ID 10) - use sqlUpdate() for the records that have an ID that is *already* in the database table (and for which you want to update the values) If you don't know in advance which ID are already in the database table, you need to read the table first (to fetch existing ID's), and use that to divide your table in 2 sets: one for already existing ID (for sqlUpdate), and one for new ID (for sqlSave) I hope that will help Cheers Julien -- View this message in context: http://r.789695.n4.nabble.com/sqlSave-tp892040p4635387.html Sent from the R help mailing list archive at Nabble.com. __ 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] Printing from R Console in colour
On 04.07.2012 08:46, amarjit chandhial wrote: Hi,I want to be able to print in colour from the R console i.e. commands (in navy) with results (in red) as on the console. From the console if I click on File - Print, the commands with results print on my printer but only in monochrome, not in colour. Similarly, if I Edit - Select All, Edit - Copy --- and paste into Word, the commands and results are only in monochrome, not in colour. How do I enable R to print commands and results in colour? Thanks You are talking about RGui in Windows? Then the answer is you cannot. Uwe Ligges [[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. __ 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] Printing from R Console in colour
I would guess that that highly depends on -what exact console you have -where exactly you paste it For example, if i copy stuff from my console (Eclipse plugin StatEt) into this Mailprogram, it is still colored. With programs similar to word, they usually have paste options that tell you what of the formatting to retain (color,font,size, etc) If you can print directly in color depends even more on the version of console your using, since, for example, mine has no such option at all. Printing in color at all depends on having a color-printer with color in it, i know it sounds dumb but check it ;) Also printer options may be set to black and white by default, so check that too. So, really, more info needed^^ On 04.07.2012, at 08:46, amarjit chandhial wrote: Hi,I want to be able to print in colour from the R console i.e. commands (in navy) with results (in red) as on the console. From the console if I click on File - Print, the commands with results print on my printer but only in monochrome, not in colour. Similarly, if I Edit - Select All, Edit - Copy --- and paste into Word, the commands and results are only in monochrome, not in colour. How do I enable R to print commands and results in colour? Thanks [[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. __ 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] 'init.win' error when installing from source
On 03.07.2012 02:13, Duncan Murdoch wrote: On 12-07-02 2:05 PM, Erin Hodgess wrote: Dear R People: I'm installing R 2-.15.1 on a Windows 32 bit machine from source. I'm getting a strange error about init.win (please see below) Does this look familiar to anyone, please? Yes, a file was missed from the tarball. If you don't want to use R-patched, you need to download RHOME/src/modules/lapack/init_win.c, e.g. from https://svn.r-project.org/R/trunk/src/modules/lapack/init_win.c There was a corrected tarball posted, but I forget the name, and can't seem to spot it. For the records, it is available from: your-CRAN-mirror/src/base/R-2/R-2.15.1-w.tar.gz Best, Uwe Ligges Duncan Murdoch Thanks, Erin Microsoft Windows [Version 6.1.7600] Copyright (c) 2009 Microsoft Corporation. All rights reserved. c:\R\R-2.15.1\src\gnuwin32make all recommended make all recommended make[1]: `MkRules' is up to date. make[4]: Nothing to be done for `svnonly'. installing C headers make[2]: Nothing to be done for `all'. make[2]: `libRblas.dll.a' is up to date. make[5]: Nothing to be done for `svnonly'. installing C headers make --no-print-directory -C ../extra/intl CFLAGS='-O3 -Wall -pedantic -mtune=core2' -f Makefile.win make --no-print-directory -C ../appl CFLAGS='-O3 -Wall -pedantic -mtune=core2' FFLAGS='-O3 -mtune=core2' -f Makefile.win make --no-print-directory -C ../nmath CFLAGS='-O3 -Wall -pedantic -mtune=core2' FFLAGS='-O3 -mtune=core2' -f Makefile.win make --no-print-directory -C ../main CFLAGS='-O3 -Wall -pedantic -mtune=core2' FFLAGS='-O3 -mtune=core2' malloc-DEFS='-DLEA_MALLOC' -f Makefile.win make --no-print-directory -C ./getline CFLAGS='-O3 -Wall -pedantic -mtune=core2' make[4]: `gl.a' is up to date. make -f Makefile.win makeMakedeps make -f Makefile.win libpcre.a make[5]: `libpcre.a' is up to date. make[4]: Nothing to be done for `all'. make -f Makefile.win makeMakedeps make -f Makefile.win libtre.a make[5]: `libtre.a' is up to date. make[4]: Nothing to be done for `all'. make[5]: `stamp' is up to date. make[5]: `liblzma.a' is up to date. make[3]: `R.dll' is up to date. cp R.dll ../../bin/i386 make[3]: Nothing to be done for `all'. make --no-print-directory -C front-ends make[2]: `COPYRIGHTS' is up to date. make --no-print-directory -C ../modules -f Makefile.win \ CFLAGS='-O3 -Wall -pedantic -mtune=core2' FFLAGS='-O3 -mtune=core2' make[5]: *** No rule to make target `init_win.o', needed by `../../../bin/i386/Rlapack.dll'. Stop. make[4]: *** [all] Error 2 make[3]: *** [all] Error 1 make[2]: *** [rmodules] Error 2 make[1]: *** [rbuild] Error 2 make: *** [all] Error 2 c:\R\R-2.15.1\src\gnuwin32 Thanks, Erin __ 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@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 contour() plot info
Using dev.hold() and dev.flush() immediately gave a huge improvement in
appearance.
Adding code to use contourLines() just once, and then plotting the saved lines
at each step gave the final polish.
Many thanks Bob Kinley
-Original Message-
From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On
Behalf Of Robert Douglas Kinley
Sent: 03 July 2012 16:17
To: Duncan Murdoch
Cc: r-help@r-project.org
Subject: Re: [R] saving contour() plot info
Many thanks for these ideas ... I'll try them, and report back
Cheers Bob Kinley
-Original Message-
From: Duncan Murdoch [mailto:murdoch.dun...@gmail.com]
Sent: 03 July 2012 15:54
To: Robert Douglas Kinley
Cc: r-help@r-project.org
Subject: Re: [R] saving contour() plot info
On 03/07/2012 10:36 AM, Robert Douglas Kinley wrote:
{ I think this message got rejected at the 1st attempt - trying again}
R 2.15.1 , windows XP
I have a very non-stationary bivariate time-series - say {xt,yt} t=1 ...
lots.
I want to do a bivariate density contour-plot of the whole series and then
step
through the series 1 second at a time plotting that second's {x,y} subset on
top of the contour
plot and losing the previous second's subset so that the effect enables you
to see
an 'animation' of how the series 'travels' through different parts of the
joint density as time passes.
The only way I've found to do this is to repeatedly call contour() before
plotting each
seconds-worth of data ... this works, but because of the time taken to do the
contour()
calculations in each step of the loop , it has an unpleasant flashing
appearance.
Generally the way to get rid of the flashing is to use dev.hold() while
plotting, then dev.flush() when done. This isn't supported on all
graphics devices.
Is it possible to run contour() just once and save the contour- plotting
info, so that in each
step of the loop I only have do the actual plotting of the contours?
It is possible to save the contour information. See ?contourLines. This
doesn't save everything (e.g. the labelling), but it might be enough for
you.
You could also try using recordPlot() after plotting the contours, then
replayPlot() when you want to reproduce them.
Duncan Murdoch
Or any other way of achieving the desired outcome.
Grateful for any guidance ... Bob Kinley
[[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.
__
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@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] Automating R script with Windows 7
I would try first without the task scheduler: Make a .bat file and run this from the command line. This way you can see what is going on without the flashing window that is opened and closed immeadiately. maviney wrote I tried to task schedule, using the following code C:\Program Files\R\R-2.12.1\bin\i386\Rscript.exe but the Rscript just flashes at me Thanks for ur assistance, its back to researching my problem -- View this message in context: http://r.789695.n4.nabble.com/Automating-R-script-with-Windows-7-tp4446693p4635380.html Sent from the R help mailing list archive at Nabble.com. __ 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] CPU usage while running R code
I am currently running an R program on a computer with 16 Gb memory (Windows7, 64 bit). When I look at the task manager, I see that only 4 out of the 8 CPUs are being used. Is this due to some missing in the R code, or should I change something to the settings of the computer? Laurent Franckx, PhD Expert VITO NV Boeretang 200, 2400 MOL, Belgium Tel. + 32 14 33 58 22 Skype: laurent.franckx laurent.fran...@vito.be Visit our website: www.vito.be/english and http://www.vito.be/transport http://www.vito.be/e-maildisclaimer __ 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] How do you impute missing data using Latent Class Model (poLCA package)
My problem is I have data with both categorial and numerical data, currently only the categorical number contains missing data, was wondering do I make a new dataframe containing only the categorical columns? How would you use Latent Class Model specifically poLCA to impute the missing data? http://www.sscnet.ucla.edu/polisci/faculty/lewis/pdf/poLCA-JSS-final.pdf The reason why I chose not to use Multiple Imputation(MI) is because according to [http://blogs.iq.harvard.edu/sss/archives/2008/09/a_handy_trick_f.shtml] MI packages assume the Multivariate Normal Distribution which may not hold for certain types of categorical and binary data. Yucel Recai, Yulei He, and Alan Zaslavsky point out in their May 2008 article in The American Statistician, naive rounding MI imputations can bias estimates, particularly when the underlying data are asymmetric or multimodal. However instead of using Yucel Recai, Yulei He, and Alan Zaslavsky's rounding strategy [ http://amstat.tandfonline.com/doi/pdf/10.1198/000313008X300912] , I opted for the Latent Class Model: http://spitswww.uvt.nl/~vermunt/ginkel2007.pdf Who are the authors of the R package poLCA. Thanks Chris [[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.
Re: [R] Please help
Dear Jun I think you should consider to submit your question to R-SIG (special interest group) about spatial data http://r-sig-geo.2731867.n2.nabble.com/ This will improve the probability you get some help. Cheers Giuseppe Calamita - Giuseppe Calamita PhD at CNR-IMAA Italian National Council of Research - Institute of Methodologies for Environmental Analysis, Tito Scalo -Potenza ITALY -- View this message in context: http://r.789695.n4.nabble.com/Please-help-tp4635370p4635382.html Sent from the R help mailing list archive at Nabble.com. __ 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 lodging package
Hi, this happen when two or more packages have functions with the same names. Using search() you can see which package are loaded in your working environment. Thus, using the name of the function will call the function of the first package you see in the search output. To use a specific function you need to specify both the package containing the function and then the function E.G. pckg:function(parameters) hope this helps - Giuseppe Calamita PhD at CNR-IMAA Italian National Council of Research - Institute of Methodologies for Environmental Analysis, Tito Scalo -Potenza ITALY -- View this message in context: http://r.789695.n4.nabble.com/problem-in-lodging-package-tp4635291p4635383.html Sent from the R help mailing list archive at Nabble.com. __ 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] problem with quilt.plot
HI to you all! I have problem with quilt.plotthis is the code: quilt.plot(long_gridded,lat_gridded,d18O_gridded,nrow=n,ncol=m,xaxt='n',yaxt='n',xlab=NA, ylab=NA,breaks=seq(d18O_gridded_1,d18O_gridded_2,length.out=number_colors),col=col,horizontal=FALSE,nlevel=number_colors - 1,legend.args=list(quote(delta^18*O~(â°)),col=black,cex=0.8,side=3,line=1)) title(xlab=quote(longitude ~ (NA^o))) axis(1,las=1,lwd=0,lwd.ticks=1.0) title(ylab=quote(latitude ~ (NA^o))) axis(2,las=2,lwd=0,lwd.ticks=1.0) Namely, when I want to execute this code it does not plot it and it does not give any error message. Therefor I don`t even know what could be wrong. Thanks for your help!!! Â Tamara [[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.
Re: [R] CPU usage while running R code
Dear Laurent, An R proces uses only one core. It is possible to use multiple cores. Have a look at he Hig-Performance and Parallel Computing task view (http://cran.freestatistics.org/web/views/HighPerformanceComputing.html) Best regards, Thierry ir. Thierry Onkelinx Instituut voor natuur- en bosonderzoek / Research Institute for Nature and Forest team Biometrie Kwaliteitszorg / team Biometrics Quality Assurance Kliniekstraat 25 1070 Anderlecht Belgium + 32 2 525 02 51 + 32 54 43 61 85 thierry.onkel...@inbo.be www.inbo.be To call in the statistician after the experiment is done may be no more than asking him to perform a post-mortem examination: he may be able to say what the experiment died of. ~ Sir Ronald Aylmer Fisher The plural of anecdote is not data. ~ Roger Brinner The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. ~ John Tukey -Oorspronkelijk bericht- Van: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] Namens Franckx Laurent Verzonden: woensdag 4 juli 2012 10:29 Aan: 'r-help@r-project.org' Onderwerp: [R] CPU usage while running R code I am currently running an R program on a computer with 16 Gb memory (Windows7, 64 bit). When I look at the task manager, I see that only 4 out of the 8 CPUs are being used. Is this due to some missing in the R code, or should I change something to the settings of the computer? Laurent Franckx, PhD Expert VITO NV Boeretang 200, 2400 MOL, Belgium Tel. + 32 14 33 58 22 Skype: laurent.franckx laurent.fran...@vito.be Visit our website: www.vito.be/english and http://www.vito.be/transport http://www.vito.be/e-maildisclaimer __ 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. * * * * * * * * * * * * * D I S C L A I M E R * * * * * * * * * * * * * Dit bericht en eventuele bijlagen geven enkel de visie van de schrijver weer en binden het INBO onder geen enkel beding, zolang dit bericht niet bevestigd is door een geldig ondertekend document. The views expressed in this message and any annex are purely those of the writer and may not be regarded as stating an official position of INBO, as long as the message is not confirmed by a duly signed 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] help with filled.contour() -
On Jul 3, 2012, at 9:01 PM, Jonathan Hughes wrote: Dear all, I can't figure out a way to have more than one plot using filled.contour() in a single plate. I tried to use layout() or par(), but the way filled.contour() is written seems to override those commands. Any suggestions would be really appreciated. Jonathan [[alternative HTML version deleted]] This is the code in filled.contour that is overriding you par() efforts. mar.orig - (par.orig - par(c(mar, las, mfrow)))$mar on.exit(par(par.orig)) w - (3 + mar.orig[2L]) * par(csi) * 2.54 layout(matrix(c(2, 1), ncol = 2L), widths = c(1, lcm(w))) par(las = las) mar - mar.orig mar[4L] - mar[2L] mar[2L] - 1 par(mar = mar) plot.new() You could rewrite the function to allocate space differently. Note this portion of the layout help page: These functions are totally incompatible with the other mechanisms for arranging plots on a device: par(mfrow),par(mfcol) and split.screen. -- David Winsemius, MD West Hartford, CT __ 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] Error in hclust?
Dear R users, I have noted a difference in the merge distances given by hclust using centroid method. For the following data: x-c(1009.9,1012.5,1011.1,1011.8,1009.3,1010.6) and using Euclidean distance, hclust using centroid method gives the following results: x.dist-dist(x) x.aah-hclust(x.dist,method=centroid) x.aah$merge [,1] [,2] [1,] -3 -6 [2,] -1 -5 [3,] -2 -4 [4,]12 [5,]34 x.aah$height [1] 0.5 0.6 0.7 0.97500 1.36875 A calculation by hand results same merges, but at different distances for latter stages: heights: 0.5 = merging 3 and 6 = G1 0.6 = merging 1 and 5 = G2 0.7 = merging 2 and 4 = G3 *1.25 = merging G1 and G2 = G4 1.92 = merging G3 and G4* It seems that hclust is not correctly computing the group centroids. Is it correct? Best regards, Mateus [[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.
[R] MS-DOS script R
Hi, I would like to run R with a script from MS-DOS. R is in My Documents and my script too. How to do? -- View this message in context: http://r.789695.n4.nabble.com/MS-DOS-script-R-tp4635398.html Sent from the R help mailing list archive at Nabble.com. __ 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] Please help
Hi Pascal,
Thanks a lot for your reply.
Then the problem become lon lat no longer be a class of dim.ncdf,
they can not be used in var.def.ncdf.
Is there any way that can made lon lat become a class of dim.ncdf?
Many thanks,
Jun
Date: Wed, 4 Jul 2012 10:39:29 +0100
From: kri...@ymail.com
Subject: Re: [R] Please help
To: chensh...@hotmail.com
CC: r-help@r-project.org
Hello,
Following lines are wrong:
x=dim.def.ncdf(Lon,degreesE,140.0251:146.6751)
y=dim.def.ncdf(Lat,degreesN,(-31.025):(-38.975))
You have 241 longitudes and 160 latitudes.
lon - seq(from=140.0251, to=146.6751, length.out=241)
lat - seq(from=-38.975, to=-31.025, length.out=160)
Regards
- Mail original -
De : Jun Chen chensh...@hotmail.com
À : r-help@r-project.org
Cc :
Envoyé le : Mercredi 4 juillet 2012 10h52
Objet : [R] Please help
Dear All,
I am a research student in environment. I have only little
programming knowledge. I am currently doing the last project about rainfall
impact on ground water quality in an area. It happens that I have to use R to
read rainfall data (3 dimension) from ASC file (*.asc), and then write them
into one NCDF file (*.nc).
I have been working very hard on study R, but I still can not fix
the problem. Could someone please as kind as point out that what the problems
are in my R script?
Firstly, this is an example of data in asc file:
NCOLS 241
NROWS 160
XLLCORNER140.00012207031
YLLCORNER -39.
CELLSIZE0.50E-01
NODATA_VALUE -99.0
166.30 160.87 155.23 149.33 143.83 138.52 133.29
128.34 123.76 119.21
115.06 110.90 107.22 103.69 100.40 97.29 94.58
92.15 90.00 87.91
86.20 84.57 83.22 81.94 81.11 80.38 79.37 78.73
79.70 79.62
---
And then this is the script I wrote:
setwd(E:/grid)
#defining dimension
x=dim.def.ncdf(Lon,degreesE,140.0251:146.6751)
y=dim.def.ncdf(Lat,degreesN,(-31.025):(-38.975))
t=dim.def.ncdf(Time,1968-01,1:12,unlim=TRUE)
#setup variable
varmr=var.def.ncdf(mr,mm,list(x,y,t),-99.00,
longname=monthly rainfall)
#create ncdf file
ncnew=create.ncdf(rainfall.nc, varmr)
#read input
files=list.files(pattern=.asc)
mrain=matrix(0:0,0,3)
for(i in files)
{rainfall=data.frame(readAsciiGrid(i))
mrain=cbind(mrain,rainfall)
}
put.var.ncdf(ncnew, mrain)
close.ncdf(ncnew)
---
[[elided Hotmail spam]]
Many thanks,
Jun
[[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.
[[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.
Re: [R] MS-DOS script R
On 04.07.2012 14:44, cindy.dol wrote: Hi, I would like to run R with a script from MS-DOS. I don't even know how to compile R for MS-DOS to run it under that OS. I am using R since 1998, but we used 64-bit Solaris or 32-bit Linux or Windows systems at the time, but no DOS left. R is in My Documents Unlikely DOS supports such a name. and my script too. How to do? Actually, I guess you are going to run R in the Windows shell. Just add the bin path of R to the PATH environment variable, go to a Windows command shell and start R via Rscript filename.R or R CMD BATCH filename.R from the folder where filename.R is located. For more details and the differences of those approaches, see the manuals. Best, Uwe Ligges -- View this message in context: http://r.789695.n4.nabble.com/MS-DOS-script-R-tp4635398.html Sent from the R help mailing list archive at Nabble.com. __ 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@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 get prediction for a variable in WinBUGS?
This came through unformatted. Please fix. I won't look into it this way.
And please send plain text only.
Best,
Uwe Ligges
On 02.07.2012 22:28, bo yu wrote:
Dear all,I am a new user of WinBUGS and need your help. After running the following
code, I got parameters of beta0 through beta4 (stats, density), but I don't know how
to get the prediction of the last value of h, the variable I set to NA and want to
model it using the following code.Does anyone can given me a hint? Any advice would
be greatly appreciated.Best Regards,Bomodel {for(i in 1: N)
{CF01[i] ~ dnorm(0, 20)CF02[i] ~ dnorm(0, 1)h[i] ~ dpois (lambda [i])log(lambda [i])
- beta0 + beta1*CF03[i] + beta2*CF02[i] + beta3*CF01[i] + beta4*IND[i]}beta0 ~
dnorm(0.0, 1.0E-6)beta1 ~ dnorm(0.0, 1.0E-6)beta2 ~ dnorm(0.0, 1.0E-6)beta3 ~
dnorm(0.0, 1.0E-6)beta4 - log(p)p ~ dunif(lower, upper)}INITSlist(beta0 = 0,
beta1 = 0, beta2 = 0, beta3 = 0, p = 0.9)DATA(LIST)list(N = 154, lower = 0.80, upper
= 0.95,h = c(1, 4, 1, 2, 1, 2, 1, 1, 1, 3, 3, 0, 0, 0, 2, 0, 1, 0, 4, 2,3, 0, 2, 1,
1, 2, 2, 2, 3, 4, 2, 3, 1, 0, 1, 3, 3, 3, 1, 0, 1,0, 5, 2, 1, 2, 1, 3, 3, 1, !
1,!
0, 2, 2, 0, 3, 0, 0, 3, 2, 2, 2,1, 0, 3, 3, 1, 1, 1, 2, 1, 0, 1, 2, 1, 2, 0,
2, 1, 0, 0, 2, 5,0, 2, 1, 0, 2, 1, 2, 2, 2, 0, 3, 2, 1, 3, 3, 3, 3, 0, 1, 3,
3,3, 1, 0, 0, 1, 2, 1, 0, 1, 4, 1, 1, 1, 1, 2, 1, 3, 0, 0, 1, 1,1, 1, 0, 2, 1,
0, 0, 1, 1, 5, 1, 1, 1, 3, 0, 1, 1, 1, 0, 2, 1,0, 3, 3, 0, 0, 1, 2, 6, NA),CF03
= c(-1.575, 0.170, -1.040, -0.010, -0.750,0.665, -0.250, 0.145, -0.345, -1.915,
-1.515,0.215, -1.040, -0.035, 0.805, -0.860, -1.775,1.725, -1.345, 1.055,
-1.935, -0.160, -0.075,-1.305, 1.175, 0.130, -1.025, -0.630, 0.065,-0.665,
0.415, -0.660, -1.145, 0.165, 0.955,-0.920, 0.250, -0.365, 0.750, 0.045,
-2.760,-0.520, -0.095, 0.700, 0.155, -0.580, -0.970,-0.685, -0.640, -0.900,
-0.250, -1.355, -1.330,0.440, -1.505, -1.715, -0.330, 1.375, -1.135,-1.285,
0.605, 0.360, 0.705, 1.380, -2.385, -1.875,-0.390, 0.770, 1.605, -0.430,
-1.120, 1.575, 0.440,-1.320, -0.540, -1.490, -1.815, -2.395, 0.305,0.735,
-0.790, -1.070, -1.085, -0.540, -0.935,-0.790, 1.400, 0.310, -1.150, -!
0.7!
25, -0.150,-0.640, 2.040, -1.180, -0.235, -0.070, -0.500,-0.750, -1.45
0, -0.235, -1.635, -0.460, -1.855,-0.925, 0.075, 2.900, -0.820, -0.170,
-0.355,-0.170, 0.595, 0.655, 0.070, 0.330, 0.395, 1.165,0.750, -0.275, -0.700,
0.880, -0.970, 1.155, 0.600,-0.075, -1.120, 1.480, -1.255, 0.255, 0.725,-1.230,
-0.760, -0.380, -0.015, -1.005, -1.605,0.435, -0.695, -1.995, 0.315, -0.385,
-0.175,-0.470, -1.215, 0.780, -1.860, -0.035, -2.700,-1.055, 1.210, 0.600,
-0.710, 0.425, 0.155, -0.525,-0.565),CF02 = c(NA, NA, NA, NA, NA, NA, NA, NA,
NA,NA, NA, NA, NA, NA, NA, 0.38, 0.06, -0.94,-0.02, -0.28, -0.78, -0.95, 2.33,
1.43, 1.24, 1.26,-0.75, -1.5, -2.09, 1.01, -0.05, 2.48, 2.48, 0.46,0.46, -0.2,
-1.11, 0.52, -0.37, 0.58, 0.86, 0.59,-0.12, -1.33, 1.4, -1.84, -1.4, -0.76,
-0.23,-1.78, -1.43, 1.2, 0.32, 1.87, 0.43, -1.71, -0.54,-1.25, -1.01, -1.98,
0.52, -1.07, -0.44, -0.24,-1.31, -2.14, -0.43, 2.47, -0.09, -1.32, -0.3,-0.99,
1.1, 0.41, 1.01, -0.19, 0.45, -0.07, -1.41,0.87, 0.68, 1.61, 0.36, -1.06,
-0.44, -0.16, 0.72,-0.69, -0.94, 0.11, 1.25, 0.33, -0.05, 0.87!
, !
-0.37,-0.2, -2.22, 0.26, -0.53, -1.59, 0.04, 0.16, -2.66,-0.21, -0.92, 0.25,
-1.36, -1.62, 0.61, -0.2, 0,1.14, 0.27, -0.64, 2.29, -0.56, -0.59, 0.44,
-0.05,0.56, 0.71, 0.32, -0.38, 0.01, -1.62, 1.74, 0.27, 0.97,1.22, -0.21,
-0.05, 1.15, 1.49, -0.15, 0.05, -0.87,-0.3, -0.08, 0.5, 0.84, -1.67, 0.69,
0.47, 0.44,-1.35, -0.24, -1.5, -1.32, -0.08, 0.76, -0.57,-0.84, -1.11, 1.94,
-0.68),CF01 = c(NA, NA, NA, NA, NA, NA, NA, NA, NA,NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA,NA, -0.117, -0.211, -0.333, -0.229, -0.272,-0.243, -0.148, 0.191,
-0.263, -0.239, -0.168,-0.381, -0.512, -0.338, -0.296, 0.067, 0.104,-0.254,
-0.167, -0.526, -0.096, -0.43, 0.013,-0.438, -0.297, -0.131, -0.098, -0.046,
-0.063,-0.194, -0.155, -0.645, -0.603, -0.374, -0.214,-0.165, -0.509, -0.171,
-0.442, -0.468, -0.289,-0.427, -0.519, -0.454, 0.046, -0.275, -0.401,-0.542,
-0.488, -0.52, -0.018, -0.551, -0.444,-0.254, -0.286, 0.048, -0.03, -0.015,
-0.219,-0.029, 0.059, 0.007, 0.157, 0.141, -0.035, 0.136,0.526, 0.113!
, 0!
.22, -0.022, -0.173, 0.021, -0.027,0.261, 0.082, -0.266, -0.284, -0.09
7, 0.097, -0.06,0.397, 0.315, 0.302, -0.026, 0.268, -0.111, 0.084,0.14, -0.073,
0.287, 0.061, 0.035, -0.022, -0.091,-0.22, -0.021, -0.17, -0.184, 0.121,
-0.192,-0.24, -0.283, -0.003, -0.45, -0.138, -0.143,0.017, -0.245, 0.003,
0.108, 0.015, -0.219, 0.09,-0.22, -0.004, -0.178, 0.396, 0.204, 0.342,
0.079,-0.034, -0.122, -0.24, -0.125, 0.382, 0.072, 0.294,0.577, 0.4, 0.213,
0.359, 0.074, 0.388, 0.253, 0.167),IND = c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1,1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1,1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
Re: [R] problem with quilt.plot
On 04.07.2012 12:39, Tamara Hunjak wrote: HI to you all! I have problem with quilt.plotthis is the code: quilt.plot(long_gridded,lat_gridded,d18O_gridded,nrow=n,ncol=m,xaxt='n',yaxt='n',xlab=NA, ylab=NA,breaks=seq(d18O_gridded_1,d18O_gridded_2,length.out=number_colors),col=col,horizontal=FALSE,nlevel=number_colors - 1,legend.args=list(quote(delta^18*O~(‰)),col=black,cex=0.8,side=3,line=1)) Strange, I get: Error: could not find function quilt.plot Or in other words: This is not reproducible at all. Please read the posting guide. Best, Uwe Ligges title(xlab=quote(longitude ~ (NA^o))) axis(1,las=1,lwd=0,lwd.ticks=1.0) title(ylab=quote(latitude ~ (NA^o))) axis(2,las=2,lwd=0,lwd.ticks=1.0) Namely, when I want to execute this code it does not plot it and it does not give any error message. Therefor I don`t even know what could be wrong. Thanks for your help!!!  Tamara [[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. __ 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] Difference between two-way ANOVA and (two-way) ANCOVA
Hi! as my subject says I am struggling with the different of a two-way ANOVA and a (two-way) ANCOVA. I found the following examples from this webpage: http://www.statmethods.net/stats/anova.html # One Way Anova (Completely Randomized Design) fit - aov(y ~ A, data=mydataframe) # Randomized Block Design (B is the blocking factor) fit - aov(y ~ A + B, data=mydataframe) # Two Way Factorial Design fit - aov(y ~ A + B + A:B, data=mydataframe) fit - aov(y ~ A*B, data=mydataframe) # same thing # Analysis of Covariance fit - aov(y ~ A + x, data=mydataframe) I) The 1. example is pretty clear. A simple on way ANOVA. II) Is it correct to say that example 2. (which is called a Randomized Block Design) is a two way ANOVA? III) Example 3 is like example 2. (in case I was right in II) ) a two way ANOVA but including an interaction term. That's why they call it here a Factorial Design. So far so good. IV) For me, the ANCOVA (last example) looks like a two-way ANOVA. So in what way is the variable x different to variable B so that it is called an ANCOVA and not an ANOVA??? I presume that from the type of data R knows whether to perform an ANCOVA or an ANOVA. V) Is it right to say that the ANCOVA example is a two-way ANCOVA? Or can a one-way ANCOVA actually exists? You see I am a bit confused especially how R distinguishes between the ANCOVA and the two-way ANOVA? I hope to find some useful answers here. Cheers! -- View this message in context: http://r.789695.n4.nabble.com/Difference-between-two-way-ANOVA-and-two-way-ANCOVA-tp4635403.html Sent from the R help mailing list archive at Nabble.com. __ 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] how to get list of files within a particular local file folder
Dear List, Say I can use getwd() and setwd() to change my working directory. How can I read in all the files within that directory using command line (like a ls() but for the path specified) Regards Ajay Websites- Technology http://decisionstats.com On Wed, Jul 4, 2012 at 3:30 PM, r-help-requ...@r-project.org wrote: r-help@r-project.org [[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.
Re: [R] how to get list of files within a particular local file folder
Dear List, Say I can use getwd() and setwd() to change my working directory. How can I read in all the files within that directory using command line (like a ls() but for the path specified) Regards Ajay Websites- Technology http://decisionstats.com [[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.
Re: [R] CPU usage while running R code
R does not parallelize its operation automatically... you have to use R code that splits the work you give it into multiple tasks. See ?parallel --- Jeff NewmillerThe . . Go Live... DCN:jdnew...@dcn.davis.ca.usBasics: ##.#. ##.#. Live Go... Live: OO#.. Dead: OO#.. Playing Research Engineer (Solar/BatteriesO.O#. #.O#. with /Software/Embedded Controllers) .OO#. .OO#. rocks...1k --- Sent from my phone. Please excuse my brevity. Franckx Laurent laurent.fran...@vito.be wrote: I am currently running an R program on a computer with 16 Gb memory (Windows7, 64 bit). When I look at the task manager, I see that only 4 out of the 8 CPUs are being used. Is this due to some missing in the R code, or should I change something to the settings of the computer? Laurent Franckx, PhD Expert VITO NV Boeretang 200, 2400 MOL, Belgium Tel. + 32 14 33 58 22 Skype: laurent.franckx laurent.fran...@vito.be Visit our website: www.vito.be/english and http://www.vito.be/transport http://www.vito.be/e-maildisclaimer __ 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@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 do you impute missing data using Latent Class Model (poLCA package)
John Kane Kingston ON Canada Who are the authors of the R package poLCA. From the R CRAN site Author: Drew Linzer, Jeffrey Lewis. Maintainer: Drew Linzer dlinzer at emory.edu FREE 3D MARINE AQUARIUM SCREENSAVER - Watch dolphins, sharks orcas on your desktop! __ 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 get list of files within a particular local file folder
How you read all ... files is up to you, as that depends both on the type of data contained in the files and on how you plan to use the data. Most likely the solution will involve using the files.list function and either lapply or a for loop, and the data will end up stored in a list of data objects. (If the data are tabular you might use the sapply function to create a single data table.) --- Jeff NewmillerThe . . Go Live... DCN:jdnew...@dcn.davis.ca.usBasics: ##.#. ##.#. Live Go... Live: OO#.. Dead: OO#.. Playing Research Engineer (Solar/BatteriesO.O#. #.O#. with /Software/Embedded Controllers) .OO#. .OO#. rocks...1k --- Sent from my phone. Please excuse my brevity. Ajay Ohri ohri2...@gmail.com wrote: Dear List, Say I can use getwd() and setwd() to change my working directory. How can I read in all the files within that directory using command line (like a ls() but for the path specified) Regards Ajay Websites- Technology http://decisionstats.com [[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. __ 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] Difference between two-way ANOVA and (two-way) ANCOVA
On Jul 4, 2012, at 15:20 , syrvn wrote: Hi! as my subject says I am struggling with the different of a two-way ANOVA and a (two-way) ANCOVA. I found the following examples from this webpage: http://www.statmethods.net/stats/anova.html # One Way Anova (Completely Randomized Design) fit - aov(y ~ A, data=mydataframe) # Randomized Block Design (B is the blocking factor) fit - aov(y ~ A + B, data=mydataframe) # Two Way Factorial Design fit - aov(y ~ A + B + A:B, data=mydataframe) fit - aov(y ~ A*B, data=mydataframe) # same thing # Analysis of Covariance fit - aov(y ~ A + x, data=mydataframe) I) The 1. example is pretty clear. A simple on way ANOVA. II) Is it correct to say that example 2. (which is called a Randomized Block Design) is a two way ANOVA? III) Example 3 is like example 2. (in case I was right in II) ) a two way ANOVA but including an interaction term. That's why they call it here a Factorial Design. So far so good. IV) For me, the ANCOVA (last example) looks like a two-way ANOVA. So in what way is the variable x different to variable B so that it is called an ANCOVA and not an ANOVA??? I presume that from the type of data R knows whether to perform an ANCOVA or an ANOVA. V) Is it right to say that the ANCOVA example is a two-way ANCOVA? Or can a one-way ANCOVA actually exists? You see I am a bit confused especially how R distinguishes between the ANCOVA and the two-way ANOVA? I hope to find some useful answers here. Well, it's not really about R, is it? Anyways, I'd call y~A+x a ONE-way ANCOVA, because it deals with the covariation of two variables (y and x) in a one-way layout. In the traditional applications, x is often independent of A (pre-randomization measurement like soil quality, etc.) so that the group means of y can be estimated as the value of the regression at the grand mean of x (adjusted means), and the mean difference between two groups is the vertical difference between the parallel regression lines. Cheers! -- View this message in context: http://r.789695.n4.nabble.com/Difference-between-two-way-ANOVA-and-two-way-ANCOVA-tp4635403.html Sent from the R help mailing list archive at Nabble.com. __ 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. -- Peter Dalgaard, Professor Center for Statistics, Copenhagen Business School Solbjerg Plads 3, 2000 Frederiksberg, Denmark Phone: (+45)38153501 Email: pd@cbs.dk Priv: pda...@gmail.com __ 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] nls problem
By interpreting the code line by line and looking at the output of the lines, I
got the
following result. It looks like it needs the fifu converted to an expression,
then
evaluated. This suggests a workaround, but doesn't answer the underlying
question about
whether this is supposed to work this way.
JN
str(fifu)
language exp(-k * x)
fifu2-as.expression(fifu)
fit2b - nls(y ~ fifu2, data = data, start = c(k = 1))
Error in lhs - rhs : non-numeric argument to binary operator
fifu2
expression(exp(-k * x))
fit2be - nls(y ~ eval(fifu2), data = data, start = c(k = 1))
fit2be
Nonlinear regression model
model: y ~ eval(fifu2)
data: data
k
1
residual sum-of-squares: 1.604e-06
Number of iterations to convergence: 1
Achieved convergence tolerance: 1.603e-06
--
Message: 13
Date: Tue, 03 Jul 2012 13:54:11 +0200
From: joerg van den hoff j.van_den_h...@hzdr.de
To: r-help@r-project.org
Subject: [R] nls problem
Message-ID: op.wgvcolzs24o...@marco.fz-rossendorf.de
Content-Type: text/plain; charset=iso-8859-15; format=flowed;
delsp=yes
hi list,
used versions: 2.12.1 and 2.14.0 under ubuntu and macosx.
I recently stumbled over a problem with `nls', which occurs if the model
is not specified explicitly but via an evaluation of a 'call' object.
simple example:
8--
nlsProblem - function (len = 5) {
#===
# purpose: to demonstrate an apparent problem with `nls' which occurs,
# if the model is specified by passing th lhs as an evaled 'call'
# object. The problem is related to the way `nls' tries to compute
# its internal variable `varIndex' which rests on the assumption that
# the dependent (y) and, possibly, the independent (x) variable
# are identified by having a length equal to the `nls' variable
# `respLength'. the problem arises when there are `varNames'
# components accidentally having this length, too.
# in the present example, setting the number of data points to
# len=2 triggers the error since the `call' object `fifu' has this
# length, too and `nls' constructs an erroneous `mf$formula' internally.
#===
#generate some data
x - seq(0, 4, len = len)
y - exp(-x)
y - rnorm(y, y, .001*y)
data - list(x = x, y = y)
#define suitable model
model - y ~ exp(-k*x)
fifu - model[[3]]
This is where my output above should be placed. JN
#this fit is fine:
fit1 - nls(model, data = data, start = c(k = 1))
print(summary(fit1))
#this fit crashes `nls' if len = 2:
fit2 - nls(y ~ eval(fifu), data = data, start = c(k = 1))
print(summary(fit2))
}
8--
to see the problem call `nlsProblem(2)'.
as explained in the above comments in the example function, I tracked it
down to the way
`nls' identifies x and y in the model expression. the problem surfaces in
the line
varIndex - n%%respLength == 0
(line 70 in the function listing from within R) which, in the case of
`fit2' in the above
example always returns a single TRUE index as long as `len != 2' (which
seems fine for the
further processing) but returns a TRUE value for the index of `fifu' as
well if `len == 2'.
question1: I'm rather sure it worked about 6 months ago with (ca.) 2.11.x
under ubuntu. have there been changes in this area?
question2: is something like the `fit2' line in the example expected to
work or not?
qeustion3: if it is not expected to work, should not the manpage include a
corresponding caveat?
question4: is there a a substitute/workaround for the `fit2' line which
still allows to specify the rhs of the model via a variable instead of a
constant (explicit) expression or function call?
the above example is of course construed but in my use case I actually
need this sort of thing. is there any chance that
the way `nls' analyzes its `model' argument can be changed to parse the
`eval(fifu)' construct correctly in all cases?
since I'm currently not subscribed to the list I'd appreciate if responses
could be Cc'ed to me directly.
thanks in advance,
joerg
__
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 get list of files within a particular local file folder
On Jul 4, 2012, at 15:24 , Ajay Ohri wrote: Dear List, Say I can use getwd() and setwd() to change my working directory. How can I read in all the files within that directory using command line (like a ls() but for the path specified) Something like all - lapply(list.files(), read.table) ? -pd Regards Ajay Websites- Technology http://decisionstats.com [[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. -- Peter Dalgaard, Professor Center for Statistics, Copenhagen Business School Solbjerg Plads 3, 2000 Frederiksberg, Denmark Phone: (+45)38153501 Email: pd@cbs.dk Priv: pda...@gmail.com __ 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] Please help
Hi Pascal,
I think I figure it out.
By doing the following, I can made lon lat become a class of
dim.ncdf:
lon - seq(from=140.0251, to=146.6751, length.out=241)
lat - seq(from=-38.975, to=-31.025, length.out=160)
x=dim.def.ncdf(Lon,degreesE,as.double(lon))
y=dim.def.ncdf(Lat,degreesN,as.double(lat))
However, after running the script, there is a error I really unable to
understand:
error at data.frame(..., check.names = FALSE) :
parameter value mean different rows 0, 31981
If you can understand what the problem is, please as kind as offer me a
direction.
Many thanks,
Jun
Date: Wed, 4 Jul 2012 10:39:29 +0100
From: kri...@ymail.com
Subject: Re: [R] Please help
To: chensh...@hotmail.com
CC: r-help@r-project.org
Hello,
Following lines are wrong:
x=dim.def.ncdf(Lon,degreesE,140.0251:146.6751)
y=dim.def.ncdf(Lat,degreesN,(-31.025):(-38.975))
You have 241 longitudes and 160 latitudes.
lon - seq(from=140.0251, to=146.6751, length.out=241)
lat - seq(from=-38.975, to=-31.025, length.out=160)
Regards
- Mail original -
De : Jun Chen chensh...@hotmail.com
À : r-help@r-project.org
Cc :
Envoyé le : Mercredi 4 juillet 2012 10h52
Objet : [R] Please help
Dear All,
I am a research student in environment. I have only little
programming knowledge. I am currently doing the last project about rainfall
impact on ground water quality in an area. It happens that I have to use R to
read rainfall data (3 dimension) from ASC file (*.asc), and then write them
into one NCDF file (*.nc).
I have been working very hard on study R, but I still can not fix
the problem. Could someone please as kind as point out that what the problems
are in my R script?
Firstly, this is an example of data in asc file:
NCOLS 241
NROWS 160
XLLCORNER140.00012207031
YLLCORNER -39.
CELLSIZE0.50E-01
NODATA_VALUE -99.0
166.30 160.87 155.23 149.33 143.83 138.52 133.29
128.34 123.76 119.21
115.06 110.90 107.22 103.69 100.40 97.29 94.58
92.15 90.00 87.91
86.20 84.57 83.22 81.94 81.11 80.38 79.37 78.73
79.70 79.62
---
And then this is the script I wrote:
setwd(E:/grid)
#defining dimension
x=dim.def.ncdf(Lon,degreesE,140.0251:146.6751)
y=dim.def.ncdf(Lat,degreesN,(-31.025):(-38.975))
t=dim.def.ncdf(Time,1968-01,1:12,unlim=TRUE)
#setup variable
varmr=var.def.ncdf(mr,mm,list(x,y,t),-99.00,
longname=monthly rainfall)
#create ncdf file
ncnew=create.ncdf(rainfall.nc, varmr)
#read input
files=list.files(pattern=.asc)
mrain=matrix(0:0,0,3)
for(i in files)
{rainfall=data.frame(readAsciiGrid(i))
mrain=cbind(mrain,rainfall)
}
put.var.ncdf(ncnew, mrain)
close.ncdf(ncnew)
---
[[elided Hotmail spam]]
Many thanks,
Jun
[[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.
[[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.
Re: [R] how to get list of files within a particular local file folder
Either dir() or list.files() Michael On Jul 4, 2012, at 8:24 AM, Ajay Ohri ohri2...@gmail.com wrote: Dear List, Say I can use getwd() and setwd() to change my working directory. How can I read in all the files within that directory using command line (like a ls() but for the path specified) Regards Ajay Websites- Technology http://decisionstats.com [[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. __ 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 check convergence of arima model
Also look at auto.arima in the forecast package. Michael On Jul 4, 2012, at 4:38 AM, Rui Barradas ruipbarra...@sapo.pt wrote: Hello, Put the fitted models in a list and then use lapply with AIC(). Like this set.seed(1) x - 1:100 + sqrt(1:100 + runif(100)) + rnorm(100) models - list() models[[1]] - arima(diff(x), order = c(1, 0, 0)) # Just to show models[[2]] - arima(diff(x), order = c(1, 0, 1)) # several models[[3]] - arima(diff(x), order = c(2, 0, 2)) # models models[[4]] - arima(diff(x), order = c(2, 0, 3)) lapply(models, AIC) Or run each model at a time through AIC(), whichever suits better. Hope this helps Rui Barradas Em 04-07-2012 10:22, Sajeeka Nanayakkara escreveu: Hi, I need to predict exchange rates using time series. So according to ACF and PACF knowledge, mixed model (ARIMA) be the appropriate and now, I need to find order of the model (p,d,q). So, several models were fitted to select the suitable model using arima(). Could you please tell me the procedure of selecting the correct order as I don't have enough time to search? Thank you. Sajeeka Nanayakkara *From:* Rui Barradas ruipbarra...@sapo.pt *To:* Sajeeka Nanayakkara nsaje...@yahoo.com *Cc:* r-help@r-project.org *Sent:* Wednesday, July 4, 2012 2:58 PM *Subject:* Re: [R] how to check convergence of arima model Hello, Inline. Em 04-07-2012 09:35, Sajeeka Nanayakkara escreveu: Hi, Sorry, since I didn't see the earlier message I resend it. I read the help page that you mentioned. But the problem is for all models, code is zero. According to that, all models were converged. Considering AIC value the best model is selected. No, arima() does not select models by AIC. That is the default behavior of ar(); arima() does NOT select models, it selects, using optim, values for the parameters of a specified model. You must choose the orders yourself. Hope this helps, Rui Barradas Is that correct procedure? The R code which was used is; model1-arima(rates,c(1,1,1)) model1 model1$code [1] 0 Sajeeka Nanayakkara *From:* Rui Barradas ruipbarra...@sapo.pt mailto:ruipbarra...@sapo.pt *To:* Sajeeka Nanayakkara nsaje...@yahoo.com mailto:nsaje...@yahoo.com *Cc:* r-help@r-project.org mailto:r-help@r-project.org *Sent:* Wednesday, July 4, 2012 2:01 PM *Subject:* Re: [R] how to check convergence of arima model Hello, Sorry, but do you read the answers to your posts? Inline. Em 04-07-2012 08:02, Sajeeka Nanayakkara escreveu: I have already fitted several models using R code; arima(rates,c(p,d,q)) And I have already answered to a question starting like this yesterday. In the mean time, the subject line changed. As I heard, best model produce the smallest AIC value, but maximum likelihood estimation procedure optimizer should converge. How to check whether maximum likelihood estimation procedure optimizer has converged or not? Yes, it was this question, the subject line was 'question'... ... And the answer was: read the manual, that I quoted, by the way. It now changed to: read the manual, period. Rui Barradas [[alternative HTML version deleted]] __ R-help@r-project.org mailto:R-help@r-project.org mailto:R-help@r-project.org mailto: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@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@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] Size of subsample in ecodist mantel()
pboot is the proportion of the sample to select, and so if it's greater than or equal to 1 you're using the entire sample rather than a subsample, so of course the limits are equal. If you look at the usage line in the help for mantel, where defaults are given, the default for pboot = 0.9. The details are discussed in the associated JSS paper that is given by citation(ecodist) Sarah On Thu, Jun 28, 2012 at 9:28 PM, nevil amos nevil.a...@gmail.com wrote: Thanks Sarah, It is not clear to me exactly how I set this value. if I enter a value for pboot then the ulim==llim X-dist(1:100) Y-dist(1:100+50*rnorm(100)) length(X) [1] 4950 print(mantel(X~Y,nperm=1000,nboot=1000,pboot=10)) mantelr pval1 pval2 pval3 llim.2.5% ulim.97.5% 0.1396906 0.001 1.000 0.001 0.1396906 0.1396906 if I do not set a value for pboot then as expected ulimllim print(mantel(X~Y,nperm=1000,nboot=1000)) mantelr pval1 pval2 pval3 llim.2.5% ulim.97.5% 0.1396906 0.001 1.000 0.001 0.1005011 0.1805416 What is the default for pboot? this does not appear to be dealt with in help for mantel On Thu, Jun 28, 2012 at 8:35 PM, Sarah Goslee sarah.gos...@gmail.com wrote: You can set it using the pboot argument. Sarah On Thursday, June 28, 2012, nevil amos wrote: What is the size of the boostrapped subsample in ecodist mantel() thanks -- Sarah Goslee http://www.functionaldiversity.org __ 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] Binary Quadratic Opt?
On Mon, Jul 02, 2012 at 06:11:37AM -0700, khris wrote: Hi, Petr, Hi Khris: If i understand the problem correctly, you have a list of (x,y) coordinates, where some sensor is located, but you do not know, which sensor is there. The database contains data for each sensor identified in some way, but you do not know the mapping between sensor identifiers from the database and the (x,y) coordinates. Is this correct? Yes. So I modelled the problem as inexact match between 2 Graphs. Since the best package on Graphs i.e. iGraph does not have any function for Graph matching I think, the problem is close to http://en.wikipedia.org/wiki/Graph_isomorphism You have estimates of the distances between the sensors using identifiers from the database. So, you know, which pairs of sensors are close. This is one graph. The other graph is the graph of closeness between the known (x,y) coordinates. You want to find a mapping between the vertices of these two graphs, which preserves edges. Yes, I agree the problem is more into Graph theoretic domain to be more precise inexact graph matching whose generalization is the Graph Isomorphism problem. The problem is more general than Graph Isomorphism. Let me define the problem more formally. We have 2 weighted undirected graphs. In one graph I know the distance of every vertex from every other vertex whereas in another graph I know only which vertices are close to a given vertex. So I know the neighboring vertices given a vertex. So the distance matrix of other Graph is incompletely known. So the question is can I find the best alignment between the 2 graphs. Ex:- G1 is know the complete distance matrix. For G2, if there are four vertices let's say (v1, v2, v3 v4) the I know edge weight (v1,v2) and (v1,v3) but have no information of edge weight(v1,v4). Similarly I know about (v2,v3) but no information about edge weights (v2,v4) or (v3,v4). So I was thinking of not to model it as general inexact Graph matching problem for then the complexity n^4. It seems the best way to model the solution is to consider only edges with are at distance of 1 unit i.e. closest edge from every vertex and not every edge from the given vertex. This will bring down the complexity from n^4 to 6*6*n^2 assuming every vertex has atmost 6 neighboring vertex. Quadratic complexity seems manageable. Ofcourse now the solution become lot more sensitive to the errors in Graph G2. Assuming best case if I have no errors in G2 i.e. for every vertex I know correctly it's closest neighbored in the rectangular grid then optimizing distance between G1 and G2 should give me best correct alignment. This seems to be the best approach under current circumstance. As far as implementation goes I think I still have to use optimization package since there are not any readily and freely available function for inexact graph matching. Petr how do you feel about it. Appreciate your feedback. Hi. I am on vacations with only occasional access to the internet. One way to solve your problem is to formulate it in a way, in which it can be solved by some existing solver. I do not know, whether this is possible. Look at self-organizing maps (Kohonen maps). They try to map points with unknown mutual relationships to a 2 dimensional grid. However, i am not sure, whether the input to this method is suitable for your problem. Another way is to write your own program. I sent some suggestions in this direction in a previous email. If you want, we can discuss this in more detail next week, when i am back from vacations. All the best, Petr. __ 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! Please recommend good books/resources on visualizing data and understanding multivariate relations...
One basic and very good one is Cleveland, W. S. (1985). The Elements of Graphing Data. Wadsworth, Inc. John Kane Kingston ON Canada -Original Message- From: comtech@gmail.com Sent: Tue, 3 Jul 2012 18:12:00 -0500 To: r-h...@stat.math.ethz.ch Subject: [R] Help! Please recommend good books/resources on visualizing data and understanding multivariate relations... Hi all, Could you please help me? I am looking for books/pointers/resources/tutorials on visualizing complex/big data and on understanding multivariate relations in complicated data. More specifically, we have categorical variables and are interested in how to visualize the categorical data and visualize data conditioned upon categorical values. Could anybody please give me some pointers? Thanks a lot! [[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. FREE 3D MARINE AQUARIUM SCREENSAVER - Watch dolphins, sharks orcas on your desktop! __ 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] RODBC tables
Dear Sir/Madam, I am desperately in need of some help. I am trying to access tables from the oracle database and inserting them into R via a data frame and I keep getting an error saying that Error in .Call(C_RODBCFetchRows, attr(channel, handle_ptr), max, buffsize, : negative length vectors are not allowed. My connection is fine and my code for this is: check-odbcConnect(dsn=,uid=*,pwd=**) There are terms instead of the *'s but I am not sure if I should disclose them because this is work-related. I have looked all over the internet and tried hundreds of solutions but have had no luck. These are the code I tried to use to get the table i called ANYTHING into R for which the negative length vectors errors came up. sqlTables(check,schema=**) nowfetchmewillwheaton-sqlFetch(check,ANYTHING) Do the errors imply that my connection is wrong, my code is wrong or that the table is not in the correct place in oracle. Please let me know if you can help or if you can give me the email address of someone who can. Kind Regards, Lorcan Treanor [[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.
Re: [R] how to check convergence of arima model
Hi, Thanks for your explanation. But still I didn't get the answer which I need. You have mentioned the code which can obtain AIC values for all models at once. My question is; the procedure of selecting the suitable model to predict using R package. If I select the model which produces the smallest AIC and maximized log likelihood values, as the suitable model is it correct? Sajeeka Nanayakkara From: Rui Barradas ruipbarra...@sapo.pt Cc: r-help@r-project.org Sent: Wednesday, July 4, 2012 3:38 PM Subject: Re: [R] how to check convergence of arima model Hello, Put the fitted models in a list and then use lapply with AIC(). Like this set.seed(1) x - 1:100 + sqrt(1:100 + runif(100)) + rnorm(100) models - list() models[[1]] - arima(diff(x), order = c(1, 0, 0)) # Just to show models[[2]] - arima(diff(x), order = c(1, 0, 1)) # several models[[3]] - arima(diff(x), order = c(2, 0, 2)) # models models[[4]] - arima(diff(x), order = c(2, 0, 3)) lapply(models, AIC) Or run each model at a time through AIC(), whichever suits better. Hope this helps Rui Barradas Em 04-07-2012 10:22, Sajeeka Nanayakkara escreveu: Hi, I need to predict exchange rates using time series. So according to ACF and PACF knowledge, mixed model (ARIMA) be the appropriate and now, I need to find order of the model (p,d,q). So, several models were fitted to select the suitable model using arima(). Could you please tell me the procedure of selecting the correct order as I don't have enough time to search? Thank you. Sajeeka Nanayakkara *From:* Rui Barradas ruipbarra...@sapo.pt *Cc:* r-help@r-project.org *Sent:* Wednesday, July 4, 2012 2:58 PM *Subject:* Re: [R] how to check convergence of arima model Hello, Inline. Em 04-07-2012 09:35, Sajeeka Nanayakkara escreveu: Hi, Sorry, since I didn't see the earlier message I resend it. I read the help page that you mentioned. But the problem is for all models, code is zero. According to that, all models were converged. Considering AIC value the best model is selected. No, arima() does not select models by AIC. That is the default behavior of ar(); arima() does NOT select models, it selects, using optim, values for the parameters of a specified model. You must choose the orders yourself. Hope this helps, Rui Barradas Is that correct procedure? The R code which was used is; model1-arima(rates,c(1,1,1)) model1 model1$code [1] 0 Sajeeka Nanayakkara *From:* Rui Barradas ruipbarra...@sapo.pt mailto:ruipbarra...@sapo.pt *Cc:* r-help@r-project.org mailto:r-help@r-project.org *Sent:* Wednesday, July 4, 2012 2:01 PM *Subject:* Re: [R] how to check convergence of arima model Hello, Sorry, but do you read the answers to your posts? Inline. Em 04-07-2012 08:02, Sajeeka Nanayakkara escreveu: I have already fitted several models using R code; arima(rates,c(p,d,q)) And I have already answered to a question starting like this yesterday. In the mean time, the subject line changed. As I heard, best model produce the smallest AIC value, but maximum likelihood estimation procedure optimizer should converge. How to check whether maximum likelihood estimation procedure optimizer has converged or not? Yes, it was this question, the subject line was 'question'... ... And the answer was: read the manual, that I quoted, by the way. It now changed to: read the manual, period. Rui Barradas [[alternative HTML version deleted]] __ R-help@r-project.org mailto:R-help@r-project.org mailto:R-help@r-project.org mailto: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. [[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.
Re: [R] how to get list of files within a particular local file folder
Dear Ajay Ohri,
Re:
Dear List,
Say I can use getwd() and setwd() to change my working directory. How can I
read in all the files within that directory using command line (like a ls()
but for the path specified)
I use this function,
# functions
getFolder - function(pat)
{
txt=file.choose()
#cat(txt,'\n')
pos=0
fname=basename(txt)
#cat(paste(\nFilename found is: ,fname,'\n'))
foln = dirname(txt)
cat(paste(\nFolder name found is: ,foln,'\n'))
drtext=dir(foln, pattern=pat, full.names = TRUE)
#cat('\n\n\n\n')
return(drtext)
}
Best,
Franklin
Dr F. Bretschneider
Dept of Biology
Kruytgebouw W711
Padualaan 8
3584 CH Utrecht
The Netherlands
__
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] About nlminb function
Hello I want to use the nlminb function but I have the objective function like characters. I can summarize the problem using the first example in the nlminb documentation. x - rnbinom(100, mu = 10, size = 10) hdev - function(par) -sum(dnbinom(x, mu = par[1], size = par[2], log = TRUE)) nlminb(c(9, 12), objective=hdev) With the last instructions we obtain appropriate results. If I have the name of the objective function stored in the element name, how can modify the next two instructions to obtain the same above results? name - 'hdev' nlminb(c(9, 12), objective=name) Thanks for the reply. Freddy. -- View this message in context: http://r.789695.n4.nabble.com/About-nlminb-function-tp4635421.html Sent from the R help mailing list archive at Nabble.com. __ 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] loop for regression
Hi, You could also use: dat1 - read.table(text= Date Stock1 Stock2 Stock3 Market 01/01/2000 1 2 3 4 01/02/2000 5 6 7 8 01/03/2000 1 2 3 4 01/04/2000 5 6 7 8 , header=TRUE, stringsAsFactors=FALSE) Stocks-dat1[,2:4] apply(Stocks,2,function(x) lm(x~Market,data=dat1)) $Stock1 Call: lm(formula = x ~ Market, data = dat1) Coefficients: (Intercept) Market -3 1 $Stock2 Call: lm(formula = x ~ Market, data = dat1) Coefficients: (Intercept) Market -2 1 $Stock3 Call: lm(formula = x ~ Market, data = dat1) Coefficients: (Intercept) Market -1 1 A.K. - Original Message - From: Akhil dua akhil.dua...@gmail.com To: r-help@r-project.org Cc: Sent: Wednesday, July 4, 2012 1:08 AM Subject: [R] loop for regression -- Forwarded message -- From: Akhil dua akhil.dua...@gmail.com Date: Wed, Jul 4, 2012 at 10:33 AM Subject: To: r-help@r-project.org Hi everyone I have data on stock prices and market indices and I need to run a seperate regression of every stock on market so I want to write a for loop so that I wont have to write codes again and again to run the regression... my data is in the format given below Date Stock1 Stock2 Stock3 Market 01/01/2000 1 2 3 4 01/02/2000 5 6 7 8 01/03/2000 1 2 3 4 01/04/2000 5 6 7 8 So can any one help me how to write this loop [[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. __ 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 nlminb function
Dear Freddy, Thank you for the explanation and the reproducible example. You can use get() as follows: nlminb(c(9, 12), objective=get(name)) Regards, Jorge.- On Wed, Jul 4, 2012 at 12:37 PM, Freddy Hernández wrote: Hello I want to use the nlminb function but I have the objective function like characters. I can summarize the problem using the first example in the nlminb documentation. x - rnbinom(100, mu = 10, size = 10) hdev - function(par) -sum(dnbinom(x, mu = par[1], size = par[2], log = TRUE)) nlminb(c(9, 12), objective=hdev) With the last instructions we obtain appropriate results. If I have the name of the objective function stored in the element name, how can modify the next two instructions to obtain the same above results? name - 'hdev' nlminb(c(9, 12), objective=name) Thanks for the reply. Freddy. -- View this message in context: http://r.789695.n4.nabble.com/About-nlminb-function-tp4635421.html Sent from the R help mailing list archive at Nabble.com. __ 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. [[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.
Re: [R] help with filled.contour() -
Hello, You can use .filled.contour (with initial dot) with par. I've tested it with one of the help page examples, reformulated to use the args of .filled.contour and it works. x - y - seq(-4*pi, 4*pi, len = 100) r - sqrt(outer(x^2, y^2, +)) z - cos(r^2)*exp(-r/(2*pi)) zlim - range(z, finite=TRUE) # levels - pretty(zlim, 20) col - heat.colors(20) # op - par(mfrow=c(1, 2)) plot.new() .filled.contour(x, y, z, levels = levels, col=col) plot.new() .filled.contour(x, y, z, levels = levels, col=col) par(op) Hope this helps, Rui Barradas Em 04-07-2012 09:39, Robert Douglas Kinley escreveu: Take a look at the code for filled.contour(). You'll find a line beginning .Internal(filledcontour( You can adapt this line and the lines around it to achieve what you want. Good luck Bob Kinley -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of Jonathan Hughes Sent: 04 July 2012 02:01 To: r-help@r-project.org Subject: [R] help with filled.contour() - Dear all, I can't figure out a way to have more than one plot using filled.contour() in a single plate. I tried to use layout() or par(), but the way filled.contour() is written seems to override those commands. Any suggestions would be really appreciated. Jonathan [[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. __ 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@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 nlminb function
Hello, Try makeFunction - function(x) eval( parse(text=x) ) name - 'hdev' nlminb(c(9, 12), objective = makeFunction(name)) Hope this helps, Rui Barradas Em 04-07-2012 17:37, Freddy Hernández escreveu: Hello I want to use the nlminb function but I have the objective function like characters. I can summarize the problem using the first example in the nlminb documentation. x - rnbinom(100, mu = 10, size = 10) hdev - function(par) -sum(dnbinom(x, mu = par[1], size = par[2], log = TRUE)) nlminb(c(9, 12), objective=hdev) With the last instructions we obtain appropriate results. If I have the name of the objective function stored in the element name, how can modify the next two instructions to obtain the same above results? name - 'hdev' nlminb(c(9, 12), objective=name) Thanks for the reply. Freddy. -- View this message in context: http://r.789695.n4.nabble.com/About-nlminb-function-tp4635421.html Sent from the R help mailing list archive at Nabble.com. __ 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@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] How to generate a correlated binary data set?
Hi. I am trying to generate a correlated binary data set. I've tried to use mvtBinaryEP, binarySimCLF, and bindata packages but none of them works in R version 2.15.1. Do you know any package to generate correlated binary covariates and work in R version 2.15.1, or how to generate it? Thanks, [[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.
Re: [R] Difference between two-way ANOVA and (two-way) ANCOVA
The usual terminology uses the number of ways to mean the number of factors (categorical or classification variables, with more than one degree of freedom per factor). The term covariate is used for continuous variables, with exactly one df. On Wed, Jul 4, 2012 at 9:20 AM, syrvn ment...@gmx.net wrote: Hi! as my subject says I am struggling with the different of a two-way ANOVA and a (two-way) ANCOVA. I found the following examples from this webpage: http://www.statmethods.net/stats/anova.html # One Way Anova (Completely Randomized Design) fit - aov(y ~ A, data=mydataframe) # Randomized Block Design (B is the blocking factor) fit - aov(y ~ A + B, data=mydataframe) # Two Way Factorial Design fit - aov(y ~ A + B + A:B, data=mydataframe) fit - aov(y ~ A*B, data=mydataframe) # same thing # Analysis of Covariance fit - aov(y ~ A + x, data=mydataframe) I) The 1. example is pretty clear. A simple on way ANOVA. II) Is it correct to say that example 2. (which is called a Randomized Block Design) is a two way ANOVA? III) Example 3 is like example 2. (in case I was right in II) ) a two way ANOVA but including an interaction term. That's why they call it here a Factorial Design. So far so good. IV) For me, the ANCOVA (last example) looks like a two-way ANOVA. So in what way is the variable x different to variable B so that it is called an ANCOVA and not an ANOVA??? I presume that from the type of data R knows whether to perform an ANCOVA or an ANOVA. V) Is it right to say that the ANCOVA example is a two-way ANCOVA? Or can a one-way ANCOVA actually exists? You see I am a bit confused especially how R distinguishes between the ANCOVA and the two-way ANOVA? I hope to find some useful answers here. Cheers! -- View this message in context: http://r.789695.n4.nabble.com/Difference-between-two-way-ANOVA-and-two-way-ANCOVA-tp4635403.html Sent from the R help mailing list archive at Nabble.com. __ 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. [[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.
Re: [R] how to check convergence of arima model
1. This is a statistical question. Please do not post further to this list. It is about R, and you have summarily dismissed all attempts to answer your R questions as unhelpful. So you need to look elsewhere. 2. Consult a statistician -- you use the statistical words, but do not understand what they mean. There often is NO single model that both maximizes likelihood and minimizes AIC (depending on the models one fits). -- Bert On Wed, Jul 4, 2012 at 7:03 AM, Sajeeka Nanayakkara nsaje...@yahoo.comwrote: Hi, Thanks for your explanation. But still I didn't get the answer which I need. You have mentioned the code which can obtain AIC values for all models at once. My question is; the procedure of selecting the suitable model to predict using R package. If I select the model which produces the smallest AIC and maximized log likelihood values, as the suitable model is it correct? Sajeeka Nanayakkara From: Rui Barradas ruipbarra...@sapo.pt Cc: r-help@r-project.org Sent: Wednesday, July 4, 2012 3:38 PM Subject: Re: [R] how to check convergence of arima model Hello, Put the fitted models in a list and then use lapply with AIC(). Like this set.seed(1) x - 1:100 + sqrt(1:100 + runif(100)) + rnorm(100) models - list() models[[1]] - arima(diff(x), order = c(1, 0, 0)) # Just to show models[[2]] - arima(diff(x), order = c(1, 0, 1)) # several models[[3]] - arima(diff(x), order = c(2, 0, 2)) # models models[[4]] - arima(diff(x), order = c(2, 0, 3)) lapply(models, AIC) Or run each model at a time through AIC(), whichever suits better. Hope this helps Rui Barradas Em 04-07-2012 10:22, Sajeeka Nanayakkara escreveu: Hi, I need to predict exchange rates using time series. So according to ACF and PACF knowledge, mixed model (ARIMA) be the appropriate and now, I need to find order of the model (p,d,q). So, several models were fitted to select the suitable model using arima(). Could you please tell me the procedure of selecting the correct order as I don't have enough time to search? Thank you. Sajeeka Nanayakkara *From:* Rui Barradas ruipbarra...@sapo.pt *Cc:* r-help@r-project.org *Sent:* Wednesday, July 4, 2012 2:58 PM *Subject:* Re: [R] how to check convergence of arima model Hello, Inline. Em 04-07-2012 09:35, Sajeeka Nanayakkara escreveu: Hi, Sorry, since I didn't see the earlier message I resend it. I read the help page that you mentioned. But the problem is for all models, code is zero. According to that, all models were converged. Considering AIC value the best model is selected. No, arima() does not select models by AIC. That is the default behavior of ar(); arima() does NOT select models, it selects, using optim, values for the parameters of a specified model. You must choose the orders yourself. Hope this helps, Rui Barradas Is that correct procedure? The R code which was used is; model1-arima(rates,c(1,1,1)) model1 model1$code [1] 0 Sajeeka Nanayakkara *From:* Rui Barradas ruipbarra...@sapo.pt mailto: ruipbarra...@sapo.pt *Cc:* r-help@r-project.org mailto:r-help@r-project.org *Sent:* Wednesday, July 4, 2012 2:01 PM *Subject:* Re: [R] how to check convergence of arima model Hello, Sorry, but do you read the answers to your posts? Inline. Em 04-07-2012 08:02, Sajeeka Nanayakkara escreveu: I have already fitted several models using R code; arima(rates,c(p,d,q)) And I have already answered to a question starting like this yesterday. In the mean time, the subject line changed. As I heard, best model produce the smallest AIC value, but maximum likelihood estimation procedure optimizer should converge. How to check whether maximum likelihood estimation procedure optimizer has converged or not? Yes, it was this question, the subject line was 'question'... ... And the answer was: read the manual, that I quoted, by the way. It now changed to: read the manual, period. Rui Barradas [[alternative HTML version deleted]] __ R-help@r-project.org mailto:R-help@r-project.org mailto:R-help@r-project.org mailto: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. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list
Re: [R] loop for regression
Please carefully read ?lm. As I previously told the OP, no looping/apply is necessary. The left hand side of the lm formula can be a matrix for which separate fits will be done on each column automatically. -- Bert On Wed, Jul 4, 2012 at 9:44 AM, arun smartpink...@yahoo.com wrote: Hi, You could also use: dat1 - read.table(text= Date Stock1 Stock2 Stock3Market 01/01/20001 2 34 01/02/20005 6 78 01/03/20001 2 34 01/04/20005 6 78 , header=TRUE, stringsAsFactors=FALSE) Stocks-dat1[,2:4] apply(Stocks,2,function(x) lm(x~Market,data=dat1)) $Stock1 Call: lm(formula = x ~ Market, data = dat1) Coefficients: (Intercept) Market -31 $Stock2 Call: lm(formula = x ~ Market, data = dat1) Coefficients: (Intercept) Market -21 $Stock3 Call: lm(formula = x ~ Market, data = dat1) Coefficients: (Intercept) Market -11 A.K. - Original Message - From: Akhil dua akhil.dua...@gmail.com To: r-help@r-project.org Cc: Sent: Wednesday, July 4, 2012 1:08 AM Subject: [R] loop for regression -- Forwarded message -- From: Akhil dua akhil.dua...@gmail.com Date: Wed, Jul 4, 2012 at 10:33 AM Subject: To: r-help@r-project.org Hi everyone I have data on stock prices and market indices and I need to run a seperate regression of every stock on market so I want to write a for loop so that I wont have to write codes again and again to run the regression... my data is in the format given below Date Stock1 Stock2 Stock3Market 01/01/2000 1 2 3 4 01/02/2000 5 6 7 8 01/03/2000 1 2 3 4 01/04/2000 5 6 7 8 So can any one help me how to write this loop [[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. __ 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. -- Bert Gunter Genentech Nonclinical Biostatistics Internal Contact Info: Phone: 467-7374 Website: http://pharmadevelopment.roche.com/index/pdb/pdb-functional-groups/pdb-biostatistics/pdb-ncb-home.htm [[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.
Re: [R] Printing from R Console in colour
As has been mentioned, the windows GUI will not do this for you, but here are some options. You can save or copy the transcript file and load it into an R syntax aware editor (e.g. emacs with ess and others) which will do the coloring/formatting for you and may be able to print with the coloring. You can create the transcript using the etxtStart function in the TeachingDemos package which can then be converted to a file with the coloring (future versions will aslo create rtf or pandoc markdown files to do something similar). You can use the knitr package (or relatives) to run your code and have the input vs. output nicely formatted. On Wed, Jul 4, 2012 at 12:46 AM, amarjit chandhial a.chandh...@btinternet.com wrote: Hi,I want to be able to print in colour from the R console i.e. commands (in navy) with results (in red) as on the console. From the console if I click on File - Print, the commands with results print on my printer but only in monochrome, not in colour. Similarly, if I Edit - Select All, Edit - Copy --- and paste into Word, the commands and results are only in monochrome, not in colour. How do I enable R to print commands and results in colour? Thanks [[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. -- Gregory (Greg) L. Snow Ph.D. 538...@gmail.com __ 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 generate a correlated binary data set?
A simple approach is to generate correlated normal data (mvrnorm function in MASS package is one way), then use a cut-off to convert them to binary. On Wed, Jul 4, 2012 at 12:25 PM, Soyeon Kim yunni0...@gmail.com wrote: Hi. I am trying to generate a correlated binary data set. I've tried to use mvtBinaryEP, binarySimCLF, and bindata packages but none of them works in R version 2.15.1. Do you know any package to generate correlated binary covariates and work in R version 2.15.1, or how to generate it? Thanks, [[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. -- Gregory (Greg) L. Snow Ph.D. 538...@gmail.com __ 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] loop for regression
On Wed, Jul 4, 2012 at 11:48 AM, Bert Gunter gunter.ber...@gene.com wrote: Please carefully read ?lm. As I previously told the OP, no looping/apply is necessary. The left hand side of the lm formula can be a matrix for which separate fits will be done on each column automatically. Which is a great option if the design matrix is constant, but suppose you want to predict each stock from every other stock in a dataset? This is a one liner with lapply: lapply(colnames(mtcars), function(n) lm(substitute(y ~ ., list(y = as.name(n))), data = mtcars)) granted, not the prettiest or most efficient thing on the planet (and falls apart for some reason with fastLm(), which I am still investigating work arounds to). The matrix outcome to lm() approach: lm(as.matrix(mtcars) ~ ., data = mtcars) yeilds perfect explanation by the variable of itself, as expected, which is not really useful. The OP did not give many details other than write a for loop. It is not clear what should be varying. If it is *just* the outcome, you are absolutely right, giving lm a matrix seems the most sensible route. Cheers, Josh -- Bert On Wed, Jul 4, 2012 at 9:44 AM, arun smartpink...@yahoo.com wrote: Hi, You could also use: dat1 - read.table(text= Date Stock1 Stock2 Stock3Market 01/01/20001 2 34 01/02/20005 6 78 01/03/20001 2 34 01/04/20005 6 78 , header=TRUE, stringsAsFactors=FALSE) Stocks-dat1[,2:4] apply(Stocks,2,function(x) lm(x~Market,data=dat1)) $Stock1 Call: lm(formula = x ~ Market, data = dat1) Coefficients: (Intercept) Market -31 $Stock2 Call: lm(formula = x ~ Market, data = dat1) Coefficients: (Intercept) Market -21 $Stock3 Call: lm(formula = x ~ Market, data = dat1) Coefficients: (Intercept) Market -11 A.K. - Original Message - From: Akhil dua akhil.dua...@gmail.com To: r-help@r-project.org Cc: Sent: Wednesday, July 4, 2012 1:08 AM Subject: [R] loop for regression -- Forwarded message -- From: Akhil dua akhil.dua...@gmail.com Date: Wed, Jul 4, 2012 at 10:33 AM Subject: To: r-help@r-project.org Hi everyone I have data on stock prices and market indices and I need to run a seperate regression of every stock on market so I want to write a for loop so that I wont have to write codes again and again to run the regression... my data is in the format given below Date Stock1 Stock2 Stock3Market 01/01/2000 1 2 3 4 01/02/2000 5 6 7 8 01/03/2000 1 2 3 4 01/04/2000 5 6 7 8 So can any one help me how to write this loop [[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. __ 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. -- Bert Gunter Genentech Nonclinical Biostatistics Internal Contact Info: Phone: 467-7374 Website: http://pharmadevelopment.roche.com/index/pdb/pdb-functional-groups/pdb-biostatistics/pdb-ncb-home.htm [[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. -- Joshua Wiley Ph.D. Student, Health Psychology Programmer Analyst II, Statistical Consulting Group University of California, Los Angeles https://joshuawiley.com/ __ 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 generate a correlated binary data set?
On 2012-07-04 11:25, Soyeon Kim wrote: Hi. I am trying to generate a correlated binary data set. I've tried to use mvtBinaryEP, binarySimCLF, and bindata packages but none of them works in R version 2.15.1. [...] You'll have to be more forthcoming about what you mean by works. What error is produced? I just tried the examples in mvtBinaryEP and in bindata and both work just fine for me. Perhaps you don't have the mvtnorm package installed? sessionInfo() R version 2.15.1 Patched (2012-06-27 r59661) Platform: x86_64-pc-mingw32/x64 (64-bit) locale: [1] LC_COLLATE=English_Canada.1252 LC_CTYPE=English_Canada.1252 [3] LC_MONETARY=English_Canada.1252 LC_NUMERIC=C [5] LC_TIME=English_Canada.1252 attached base packages: [1] stats graphics grDevices utils datasets methods base other attached packages: [1] bindata_0.9-18e1071_1.6 class_7.3-4 mvtBinaryEP_1.0.1 [5] mvtnorm_0.9-9992 loaded via a namespace (and not attached): [1] tools_2.15.1 Peter Ehlers __ 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 quilt.plot
On 2012-07-04 03:39, Tamara Hunjak wrote: HI to you all! I have problem with quilt.plotthis is the code: quilt.plot(long_gridded,lat_gridded,d18O_gridded,nrow=n,ncol=m,xaxt='n',yaxt='n',xlab=NA, ylab=NA,breaks=seq(d18O_gridded_1,d18O_gridded_2,length.out=number_colors),col=col,horizontal=FALSE,nlevel=number_colors - 1,legend.args=list(quote(delta^18*O~(‰)),col=black,cex=0.8,side=3,line=1)) title(xlab=quote(longitude ~ (NA^o))) axis(1,las=1,lwd=0,lwd.ticks=1.0) title(ylab=quote(latitude ~ (NA^o))) axis(2,las=2,lwd=0,lwd.ticks=1.0) Namely, when I want to execute this code it does not plot it and it does not give any error message. Therefor I don`t even know what could be wrong. 1. You should tell us that quilt.plot is in the fields package. 2. Do the examples in the documentation produce plots for you? 3. Don't include unnecessary commands like your axis labelling, unless they're relevant to the problem. Reduce the quilt.plot call to the minimum that causes the problem. 4. Try to provide _reproducible_ data to use in the command that's giving you problems. 5. I see some strange characters in the 'legend' specification - maybe that's just due to sending HTML mail which you should not do. Peter Ehlers Thanks for your help!!!  Tamara [[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.
Re: [R] How to generate a correlated binary data set?
Thanks for the response. It turns out it was not a problem of R but the problem on setting options in R studio :( Thank you! On Wed, Jul 4, 2012 at 2:25 PM, Peter Ehlers ehl...@ucalgary.ca wrote: On 2012-07-04 11:25, Soyeon Kim wrote: Hi. I am trying to generate a correlated binary data set. I've tried to use mvtBinaryEP, binarySimCLF, and bindata packages but none of them works in R version 2.15.1. [...] You'll have to be more forthcoming about what you mean by works. What error is produced? I just tried the examples in mvtBinaryEP and in bindata and both work just fine for me. Perhaps you don't have the mvtnorm package installed? sessionInfo() R version 2.15.1 Patched (2012-06-27 r59661) Platform: x86_64-pc-mingw32/x64 (64-bit) locale: [1] LC_COLLATE=English_Canada.1252 LC_CTYPE=English_Canada.1252 [3] LC_MONETARY=English_Canada.**1252 LC_NUMERIC=C [5] LC_TIME=English_Canada.1252 attached base packages: [1] stats graphics grDevices utils datasets methods base other attached packages: [1] bindata_0.9-18e1071_1.6 class_7.3-4 mvtBinaryEP_1.0.1 [5] mvtnorm_0.9-9992 loaded via a namespace (and not attached): [1] tools_2.15.1 Peter Ehlers [[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.
Re: [R] ggplot2: legend
I played around with this for a while with no success at all. I'd suggest posting the question on the ggplot2 newsgroup in Google Groups John Kane Kingston ON Canada -Original Message- From: thorn.tha...@rdls.nestle.com Sent: Tue, 3 Jul 2012 18:50:52 +0200 To: r-help@r-project.org Subject: [R] ggplot2: legend Dear all, I produced the following graph with ggplot which is almost fine, yet I don't like that the legend for Means and Observations includes a line, though no line is used in the plot for those two (the line for Overall Mean on the other hand is wanted): library(ggplot2) ddf - data.frame(x = factor(rep(LETTERS[1:2], 5)), y = rnorm(10)) p - ggplot(ddf, aes(x = x, y = y)) p + geom_point(aes(colour=Observations, shape=Observations)) + stat_summary(aes(colour = Means, shape =Means), fun.y = mean, fun.ymin = min, fun.ymax = max, geom = point) + geom_hline(aes(yintercept = mean(y), linetype = Overall Mean), show_guide = TRUE) I tried to map the linetype in geom_point (and stat_summary) to NULL and I tried to set it to blank but neither worked. Is it furthermore possible to combine the two legends? My preferred plot would have two symbols with different colours and shapes for Means and Observations (but no line) and directly below that a line for Overall Mean (that is all the three items in one legend)? For that I tried to assign the same names to the legends but this did not work either. So any help would be highly appreciated. Kind Regards, Thorn Thaler Mathematician Applied Mathematics Nestec Ltd, Nestlé Research Center PO Box 44 CH-1000 Lausanne 26 Phone: +41 21 785 8220 Fax: +41 21 785 9486 __ 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. FREE ONLINE PHOTOSHARING - Share your photos online with your friends and family! Visit http://www.inbox.com/photosharing to find out more! __ 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] loop for regression
HI Bert Thanks for the reply. You are right. In the case if it was a matrix, then lm automatically fits each column. Stocks-dat1[,2:4] is.matrix(Stocks) [1] FALSE is.data.frame(Stocks) [1] TRUE #Not working models-lm(Stocks~Market,data=dat1) Error in model.frame.default(formula = Stocks ~ Market, data = dat1, drop.unused.levels = TRUE) : #Works Stocks1-as.matrix(Stocks) models1-lm(Stocks1~Market,data=dat1) models1 Call: lm(formula = Stocks1 ~ Market, data = dat1) Coefficients: Stock1 Stock2 Stock3 (Intercept) -3 -2 -1 Market 1 1 1 A.K. From: Bert Gunter gunter.ber...@gene.com To: arun smartpink...@yahoo.com Cc: Akhil dua akhil.dua...@gmail.com; R help r-help@r-project.org Sent: Wednesday, July 4, 2012 2:48 PM Subject: Re: [R] loop for regression Please carefully read ?lm. As I previously told the OP, no looping/apply is necessary. The left hand side of the lm formula can be a matrix for which separate fits will be done on each column automatically. -- Bert On Wed, Jul 4, 2012 at 9:44 AM, arun smartpink...@yahoo.com wrote: Hi, You could also use: dat1 - read.table(text= Date Stock1 Stock2 Stock3 Market 01/01/2000 1 2 3 4 01/02/2000 5 6 7 8 01/03/2000 1 2 3 4 01/04/2000 5 6 7 8 , header=TRUE, stringsAsFactors=FALSE) Stocks-dat1[,2:4] apply(Stocks,2,function(x) lm(x~Market,data=dat1)) $Stock1 Call: lm(formula = x ~ Market, data = dat1) Coefficients: (Intercept) Market -3 1 $Stock2 Call: lm(formula = x ~ Market, data = dat1) Coefficients: (Intercept) Market -2 1 $Stock3 Call: lm(formula = x ~ Market, data = dat1) Coefficients: (Intercept) Market -1 1 A.K. - Original Message - From: Akhil dua akhil.dua...@gmail.com To: r-help@r-project.org Cc: Sent: Wednesday, July 4, 2012 1:08 AM Subject: [R] loop for regression -- Forwarded message -- From: Akhil dua akhil.dua...@gmail.com Date: Wed, Jul 4, 2012 at 10:33 AM Subject: To: r-help@r-project.org Hi everyone I have data on stock prices and market indices and I need to run a seperate regression of every stock on market so I want to write a for loop so that I wont have to write codes again and again to run the regression... my data is in the format given below Date Stock1 Stock2 Stock3 Market 01/01/2000 1 2 3 4 01/02/2000 5 6 7 8 01/03/2000 1 2 3 4 01/04/2000 5 6 7 8 So can any one help me how to write this loop [[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. __ 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. -- Bert Gunter Genentech Nonclinical Biostatistics Internal Contact Info: Phone: 467-7374 Website: http://pharmadevelopment.roche.com/index/pdb/pdb-functional-groups/pdb-biostatistics/pdb-ncb-home.htm __ 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] loop for regression
Apologies to Bert, I see now that market is a variable in the dataset, not a generic name for market indicators. lm() with a matrix as the outcome seems the best approach here. Josh On Wed, Jul 4, 2012 at 12:12 PM, Joshua Wiley jwiley.ps...@gmail.com wrote: On Wed, Jul 4, 2012 at 11:48 AM, Bert Gunter gunter.ber...@gene.com wrote: Please carefully read ?lm. As I previously told the OP, no looping/apply is necessary. The left hand side of the lm formula can be a matrix for which separate fits will be done on each column automatically. Which is a great option if the design matrix is constant, but suppose you want to predict each stock from every other stock in a dataset? This is a one liner with lapply: lapply(colnames(mtcars), function(n) lm(substitute(y ~ ., list(y = as.name(n))), data = mtcars)) granted, not the prettiest or most efficient thing on the planet (and falls apart for some reason with fastLm(), which I am still investigating work arounds to). The matrix outcome to lm() approach: lm(as.matrix(mtcars) ~ ., data = mtcars) yeilds perfect explanation by the variable of itself, as expected, which is not really useful. The OP did not give many details other than write a for loop. It is not clear what should be varying. If it is *just* the outcome, you are absolutely right, giving lm a matrix seems the most sensible route. Cheers, Josh -- Bert On Wed, Jul 4, 2012 at 9:44 AM, arun smartpink...@yahoo.com wrote: Hi, You could also use: dat1 - read.table(text= Date Stock1 Stock2 Stock3Market 01/01/20001 2 34 01/02/20005 6 78 01/03/20001 2 34 01/04/20005 6 78 , header=TRUE, stringsAsFactors=FALSE) Stocks-dat1[,2:4] apply(Stocks,2,function(x) lm(x~Market,data=dat1)) $Stock1 Call: lm(formula = x ~ Market, data = dat1) Coefficients: (Intercept) Market -31 $Stock2 Call: lm(formula = x ~ Market, data = dat1) Coefficients: (Intercept) Market -21 $Stock3 Call: lm(formula = x ~ Market, data = dat1) Coefficients: (Intercept) Market -11 A.K. - Original Message - From: Akhil dua akhil.dua...@gmail.com To: r-help@r-project.org Cc: Sent: Wednesday, July 4, 2012 1:08 AM Subject: [R] loop for regression -- Forwarded message -- From: Akhil dua akhil.dua...@gmail.com Date: Wed, Jul 4, 2012 at 10:33 AM Subject: To: r-help@r-project.org Hi everyone I have data on stock prices and market indices and I need to run a seperate regression of every stock on market so I want to write a for loop so that I wont have to write codes again and again to run the regression... my data is in the format given below Date Stock1 Stock2 Stock3Market 01/01/2000 1 2 3 4 01/02/2000 5 6 7 8 01/03/2000 1 2 3 4 01/04/2000 5 6 7 8 So can any one help me how to write this loop [[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. __ 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. -- Bert Gunter Genentech Nonclinical Biostatistics Internal Contact Info: Phone: 467-7374 Website: http://pharmadevelopment.roche.com/index/pdb/pdb-functional-groups/pdb-biostatistics/pdb-ncb-home.htm [[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. -- Joshua Wiley Ph.D. Student, Health Psychology Programmer Analyst II, Statistical Consulting Group University of California, Los Angeles https://joshuawiley.com/ -- Joshua Wiley Ph.D. Student, Health Psychology Programmer Analyst II, Statistical Consulting Group University of California, Los Angeles https://joshuawiley.com/ __ 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
[R] a problem about WLS
I was asked to do a WLS estimation, so I thought of lm() with weights like wls=lm(Y~X-1,weight=INC) however, it gives different result as below code, which use the formula of WLS y-Y*INC^-0.5 x-X*INC^-0.5 wls=lm(y~x-1) Can anybody explain to me why the first code can not give the right answer? Many thanks -- View this message in context: http://r.789695.n4.nabble.com/a-problem-about-WLS-tp4635446.html Sent from the R help mailing list archive at Nabble.com. __ 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] Please help
Hello,
A short code to include real monthly calendar to your data, similar to
the one inside NCEP2 reanalysis.
library(ncdf)
lon - seq(from=140.0251, to=146.6751, length.out=241)
lat - seq(from=-38.975, to=-31.025, length.out=160)
x=dim.def.ncdf(Lon,degreesE,as.double(lon))
y=dim.def.ncdf(Lat,degreesN,as.double(lat))
y1 = 1800# start of the period
y2 = 2012# end of the period
year - seq(y1,y2,1)
day -
c(31,28,31,30,31,30,31,31,30,31,30,31)%*%matrix(1,1,length(year));
day[2,leap.year(year)] - 29; day - as.vector(day)
hour - day*24; hour - c(hour[1],cumsum(hour[1:length(hour)-1]))
year - rep(year,each=12)
syntime - hour[year==1968] #change the year or change this line to
include more years
t=dim.def.ncdf(Time,hours since 1800-1-1 00:00:00,syntime,unlim=TRUE)
And your mrain matrix should have 3 dimensions (lon x lat x time),
before to save it as a NetCDF file.
Regards.
Le 04/07/2012 23:13, Jun Chen a écrit :
Hi Pascal,
I think I figure it out.
By doing the following, I can made lon lat become a class of
dim.ncdf:
lon - seq(from=140.0251, to=146.6751, length.out=241)
lat - seq(from=-38.975, to=-31.025, length.out=160)
x=dim.def.ncdf(Lon,degreesE,as.double(lon))
y=dim.def.ncdf(Lat,degreesN,as.double(lat))
However, after running the script, there is a error I really unable to
understand:
error at data.frame(..., check.names = FALSE) :
parameter value mean different rows 0, 31981
If you can understand what the problem is, please as kind as offer me a
direction.
Many thanks,
Jun
Date: Wed, 4 Jul 2012 10:39:29 +0100
From: kri...@ymail.com
Subject: Re: [R] Please help
To: chensh...@hotmail.com
CC: r-help@r-project.org
Hello,
Following lines are wrong:
x=dim.def.ncdf(Lon,degreesE,140.0251:146.6751)
y=dim.def.ncdf(Lat,degreesN,(-31.025):(-38.975))
You have 241 longitudes and 160 latitudes.
lon - seq(from=140.0251, to=146.6751, length.out=241)
lat - seq(from=-38.975, to=-31.025, length.out=160)
Regards
- Mail original -
De : Jun Chen chensh...@hotmail.com
À : r-help@r-project.org
Cc :
Envoyé le : Mercredi 4 juillet 2012 10h52
Objet : [R] Please help
Dear All,
I am a research student in environment. I have only little
programming knowledge. I am currently doing the last project about rainfall
impact on ground water quality in an area. It happens that I have to use R to
read rainfall data (3 dimension) from ASC file (*.asc), and then write them
into one NCDF file (*.nc).
I have been working very hard on study R, but I still can not fix
the problem. Could someone please as kind as point out that what the problems
are in my R script?
Firstly, this is an example of data in asc file:
NCOLS 241
NROWS 160
XLLCORNER140.00012207031
YLLCORNER -39.
CELLSIZE0.50E-01
NODATA_VALUE -99.0
166.30 160.87 155.23 149.33 143.83 138.52 133.29
128.34 123.76 119.21
115.06 110.90 107.22 103.69 100.40 97.29 94.58
92.15 90.00 87.91
86.20 84.57 83.22 81.94 81.11 80.38 79.37 78.73
79.70 79.62
---
And then this is the script I wrote:
setwd(E:/grid)
#defining dimension
x=dim.def.ncdf(Lon,degreesE,140.0251:146.6751)
y=dim.def.ncdf(Lat,degreesN,(-31.025):(-38.975))
t=dim.def.ncdf(Time,1968-01,1:12,unlim=TRUE)
#setup variable
varmr=var.def.ncdf(mr,mm,list(x,y,t),-99.00,
longname=monthly rainfall)
#create ncdf file
ncnew=create.ncdf(rainfall.nc, varmr)
#read input
files=list.files(pattern=.asc)
mrain=matrix(0:0,0,3)
for(i in files)
{rainfall=data.frame(readAsciiGrid(i))
mrain=cbind(mrain,rainfall)
}
put.var.ncdf(ncnew, mrain)
close.ncdf(ncnew)
---
[[elided Hotmail spam]]
Many thanks,
Jun
[[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.
__
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] EM algorithm to find MLE of coeff in mixed effects model
Dear Paul,
Thank you for your suggestion. I was moved by the fact that people are so
nice to help learners and ask for nothing.
With your help, I made some changes to my code and add more comments, try
to make things more clear.
If this R email list allow me to upload a pdf file to illustrate the
formula, it will be great.
The reason I do not use lme4 is that later I plan to use this for a more
complicated model (Y,T), Y is the same response of mixed model and T is a
drop out process.
(Frankly I did it for the more complicated model first and got NaN hessian
after one iteration, so I turned to the simple model.)
In the code, the m loop is the EM iterations. About efficiency, thank you
again for pointing it out. I compared sapply and for loop, they are tie and
for loop is even better sometimes.
In a paper by Bates, he mentioned that the asymptotic properties of beta is
kind of independent of the other two parameters. But it seems that paper
was submitted but I did not find it on google scholar.
Not sure if it is the reason for my results. That is why I hope some expert
who have done some simulation similar may be willing to give some
suggestion.
I may turn to other methods to calculate the conditional expection of the
latent variable as the same time.
Modified code is as below:
# library(PerformanceAnalytics)
# install.packages(gregmisc)
# install.packages(statmod)
library(gregmisc)
library(MASS)
library(statmod)
## function to calculate loglikelihood
loglike - function(datai = data.list[[i]], vvar = var.old,
beta = beta.old, psi = psi.old) {
temp1 - -0.5 * log(det(vvar * diag(nrow(datai$Zi)) + datai$Zi %*%
psi %*% t(datai$Zi)))
temp2 - -0.5 * t(datai$yi - datai$Xi %*% beta) %*% solve(vvar *
diag(nrow(datai$Zi)) + datai$Zi %*% psi %*% t(datai$Zi)) %*%
(datai$yi - datai$Xi %*% beta)
temp1 + temp2
}
## functions to evaluate the conditional expection, saved as Efun v2.R
## Eh1new=E(bibi')
## Eh2new=E(X'(y-Zbi))
## Eh3new=E(bi'Z'Zbi)
## Eh4new=E(Y-Xibeta-Zibi)'(Y-Xibeta-Zibi)
## Eh5new=E(Xibeta-yi)'Zibi
## Eh6new=E(bi)
Eh1new - function(datai = data.list[[i]], weights.m = weights.mat) {
bi - datai$b
tempb - bi * rep(sqrt(weights.m[, 1] * weights.m[, 2]),
2) #quadratic b, so needsqrt
t(tempb) %*% tempb/pi
} # end of Eh1
# Eh2new=E(X'(y-Zbi))
Eh2new - function(datai = data.list[[i]], weights.m = weights.mat) {
bi - datai$b
tempb - bi * rep(weights.m[, 1] * weights.m[, 2], 2)
tt - t(datai$Xi) %*% datai$yi - t(datai$Xi) %*% datai$Zi %*%
colSums(tempb)/pi
c(tt)
} # end of Eh2
# Eh3new=E(b'Z'Zbi)
Eh3new - function(datai = data.list[[i]], weights.m = weights.mat) {
bi - datai$b
tempb - bi * rep(sqrt(weights.m[, 1] * weights.m[, 2]),
2) #quadratic b, so need sqrt
sum(sapply(as.list(data.frame(datai$Zi %*% t(tempb))),
function(s) {
sum(s^2)
}))/pi
} # end of Eh3
# Eh4new=E(Y-Xibeta-Zibi)'(Y-Xibeta-Zibi)
Eh4new - function(datai = data.list[[i]], weights.m = weights.mat,
beta = beta.old) {
bi - datai$b
tt - sapply(as.list(ata.frame(c(datai$yi - datai$Xi %*%
beta) - datai$Zi %*% t(bi))), function(s) {
sum(s^2)
}) * (weights.m[, 1] * weights.m[, 2])
sum(tt)/pi
} # end of Eh4
Eh4newv2 - function(datai = data.list[[i]], weights.m = weights.mat,
beta = beta.old) {
bi - datai$b
v - weights.m[, 1] * weights.m[, 2]
temp - c()
for (i in 1:length(v)) {
temp[i] - sum(((datai$yi - datai$Xi %*% beta - datai$Zi %*%
bi[i, ]) * sqrt(v[i]))^2)
}
sum(temp)/pi
} # end of Eh4
# Eh5new=E(Xibeta-yi)'Zib
Eh5new - function(datai = data.list[[i]], weights.m = weights.mat,
beta = beta.old) {
bi - datai$b
tempb - bi * rep(weights.m[, 1] * weights.m[, 2], 2)
t(datai$Xi %*% beta - datai$yi) %*% (datai$Zi %*% t(tempb))
sum(t(datai$Xi %*% beta - datai$yi) %*% (datai$Zi %*%
t(tempb)))/pi
}
Eh6new - function(datai = data.list[[i]], weights.m = weights.mat) {
bi - datai$b
tempw - weights.m[, 1] * weights.m[, 2]
for (i in 1:nrow(bi)) {
bi[i, ] - bi[i, ] * tempw[i]
}
colMeans(bi)/pi
} # end of Eh6
## the main R script
### initial the values and generate the data
##
n - 100 #number of
observations
beta - c(-0.5, 1) #true coefficient of fixed
effects
vvar - 2 #sigma^2=2 #true error variance
of epsilon
psi - matrix(c(1, 0.2, 0.2, 1), 2, 2) #covariance matrix
of random effects b
b.m.true - mvrnorm(n = n, mu = c(0, 0), Sigma = psi) #100*2 matrix, each
row is the b_i
Xi - cbind(rnorm(7, mean = 2, sd = 0.5), log(2:8))
Zi - Xi
y.m - matrix(NA, nrow = n, ncol = nrow(Xi))#100*7, each row is
a y_i Zi=Xi
for( i in 1:n)
{
y.m[i,]=mvrnorm(1,mu=(Xi%*%beta+Zi%*%b.m.true[i,]),vvar*diag(nrow(Xi)))
}
b.list -
Re: [R] a problem about WLS
On Thu, Jul 5, 2012 at 11:40 AM, jacquesliu jacques...@gmail.com wrote: I was asked to do a WLS estimation, so I thought of lm() with weights like wls=lm(Y~X-1,weight=INC) however, it gives different result as below code, which use the formula of WLS y-Y*INC^-0.5 x-X*INC^-0.5 wls=lm(y~x-1) Can anybody explain to me why the first code can not give the right answer? The two examples are using the opposite weights. To replicate the answer from the second approach with lm, use weight=INC^-1 In lm(), observations with high weights have more, um, weight. They influence the results more than observations with low weight. Your code does the opposite. -thomas -- Thomas Lumley Professor of Biostatistics University of Auckland __ 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] Return
Hello Every one I have data on Stock prices and I want to calculate the return on all the stocks and then replace all the stock prices with the returns can any one tell me how to do My data is in the format given below Date Stock1 Stock2 Stock3 01/01/20001 2 3 01/02/20005 6 7 01/03/20001 2 3 01/04/20005 6 7 Thanks [[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.
Re: [R] Return
Perhaps something along the lines of (untested) t(apply(yourdata[, -1], 2, function(x) diff(log(x is what you are looking for. See ?apply, ?diff and ?t for more information. HTH, Jorge.- On Thu, Jul 5, 2012 at 12:58 AM, Akhil dua wrote: Hello Every one I have data on Stock prices and I want to calculate the return on all the stocks and then replace all the stock prices with the returns can any one tell me how to do My data is in the format given below Date Stock1 Stock2 Stock3 01/01/20001 2 3 01/02/20005 6 7 01/03/20001 2 3 01/04/20005 6 7 Thanks [[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. [[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.
