Re: [R] Input and output of time series data - any function or packages that helps?
On Oct 1, 2012, at 4:38 AM, jpm miao miao...@gmail.com wrote: Hello, I work with time series data. From time to time I run programs to produce results that are in time series form (e.g., quarterly or monthly data). After a few days I might need to access part of the results and to run another program. Is there any function or package (like dataframe or zoo?) that might help so that I don't need to copy the results manually to a csv or xls file? If the data are not time series (just indexed by 1, 2,3), is there any function that can help? Many --start with write.csv() Michael Thanks, Miao [[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] Input and output of time series data - any function or packages that helps?
I think this is a very vague question. If you are exporting to some other software, csv is probably a very good choice, unless there is something more specific about that software. As for manual, I think we left the scribes in the middle ages. R is eminently scriptable for automating tasks, but I know nothing of your workflow, so I think the ball is in your court. --- 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. jpm miao miao...@gmail.com wrote: Hello, I work with time series data. From time to time I run programs to produce results that are in time series form (e.g., quarterly or monthly data). After a few days I might need to access part of the results and to run another program. Is there any function or package (like dataframe or zoo?) that might help so that I don't need to copy the results manually to a csv or xls file? If the data are not time series (just indexed by 1, 2,3), is there any function that can help? Thanks, Miao [[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] merge.zoo returns unmatched dates
Sorry for the lack of reproducible data, but this seems to be a problem inherent to my dataset and I can't figure out where the issue is. I have several data frames set up as a time series with identical POSIXct date formats. If I keep the original data in data frame format and merge them using base merge- everything is perfect and everyone is happy. If I transform the data frames to zoo objects, and then do a merge.zoo- the data seem to become uncoupled from the original data. Even more unusual is that some dates in the new merged data set are prior to the original data set. I've attempted bellow to show what this looks like, and I hope someone has a suggestion as to what may be causing the problem. Here is one set of data in data.frame format head(Vup) Date Velocity_m/s 1 2010-01-21 07:42:00 1.217943 2 2010-01-21 07:43:00 1.624395 3 2010-01-21 07:44:00 1.526379 4 2010-01-21 07:45:00 1.456831 5 2010-01-21 07:46:00 1.245390 6 2010-01-21 07:47:00 1.374330 str(Vup) 'data.frame':7168 obs. of 2 variables: $ Date: POSIXct, format: 2010-01-21 07:42:00 2010-01-21 07:43:00 ... $ Velocity_m/s: num 1.22 1.62 1.53 1.46 1.25 ... And here is a second in data.frame format: head(PAS) Date PAS 1 2010-01-21 05:01:00 0.0013938 2 2010-01-21 05:02:00 0.0015331 3 2010-01-21 05:03:00 0.0016725 4 2010-01-21 05:04:00 0.0016725 5 2010-01-21 05:05:00 0.0012265 6 2010-01-21 05:06:00 0.0015889 str(PAS) 'data.frame':5520 obs. of 2 variables: $ Date : POSIXct, format: 2010-01-21 05:01:00 2010-01-21 05:02:00 ... $ PAS: num 0.00139 0.00153 0.00167 0.00167 0.00123 ... Using zoo: PASmin-zoo(as.matrix(PAS[,2]),as.POSIXct(PAS[,1],format=%Y-%m-%d %H:%M:%S,tz=UTC)) str(PASmin) zoo series from 2010-01-21 05:01:00 to 2010-01-27 13:01:00 Data: num [1:5520, 1] 0.00139 0.00153 0.00167 0.00167 0.00123 ... - attr(*, dimnames)=List of 2 ..$ : NULL ..$ : chr PAS Index: POSIXct[1:5520], format: 2010-01-21 05:01:00 2010-01-21 05:02:00 2010-01-21 05:03:00 ... ADP_UPmin-zoo(as.matrix(Vup[,2]),as.POSIXct(Vup[,1], format=%Y-%m-%d %H:%M,tz=UTC)) str(ADP_UPmin) zoo series from 2010-01-21 07:42:00 to 2010-01-26 20:12:00 Data: num [1:7168, 1] 1.22 1.62 1.53 1.46 1.25 ... - attr(*, dimnames)=List of 2 ..$ : NULL ..$ : chr UP_Velocity_m/s Index: POSIXct[1:7168], format: 2010-01-21 07:42:00 2010-01-21 07:43:00 2010-01-21 07:44:00 ... And if I merge the two zoo objects I get this: M-merge(ADP_UPmin,PASmin) head(M) UP_Velocity_m/s PAS 2010-01-20 21:01:00 NA 0.0013938 2010-01-20 21:02:00 NA 0.0015331 2010-01-20 21:03:00 NA 0.0016725 2010-01-20 21:04:00 NA 0.0016725 2010-01-20 21:05:00 NA 0.0012265 2010-01-20 21:06:00 NA 0.0015889 zoo series from 2010-01-20 21:01:00 to 2010-01-27 05:01:00 Data: num [1:8499, 1:2] NA NA NA NA NA NA NA NA NA NA ... - attr(*, dimnames)=List of 2 ..$ : NULL ..$ : chr [1:2] UP_Velocity_m/s PAR Index: POSIXct[1:8499], format: 2010-01-20 21:01:00 2010-01-20 21:02:00 2010-01-20 21:03:00 ... For some reason I can not figure out, even though both the PAS data frame and PAS zoo object starts at 2010-01-21 05:01:00, once merged the PAS data starts a day earlier at 2010-01-20 21:01:00. The actual numeric data looks good, but both variables have no come uncoupled from the time series dates (The Velocity data is similarity uncoupled). And as stated before, doing an non-zoo merge on the data.frame data works fine. Anyone got any ideas what's going on? [[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] Error in JohnsonFit
Dear Sir/Madam, As mentioned in help, R reports error while fitting Johnson distribution by method of moments. I have used moments from Weibull distribution and hence it is well within the feasible area and try to fit Johnson distribution using method of moments. As for example, shape=0.5 # Shape Parameter of IID Weibull Stress # scale=1 # Scale Parameter of IID Weibull Stress # # kth order raw moment of Weibull # k_raw= function(k) (scale^k)*(gamma(1+k/shape)) # Cumulants # k1=k_raw(1) k2=k_raw(2)-(k_raw(1))^2 k3=k_raw(3)-3*k_raw(2)*k_raw(1)+3*(k_raw(1))^3 -(k_raw(1))^3 k4=k_raw(4)-4*k3*k1-3*(k2)^2 - 6*(k2)*(k1)^2-(k1)^4 # Central Moments # theta=c( k1, k2 , k3 , k4+3*(k2)^2) parms-JohnsonFit(theta,moment=use) Error in JohnsonFit(theta, moment = use) : Couldn't do an Sb fit I will be thankful if you please let me know elaborately the exact difficulty which causes such errors. Regards, Prajamitra Bhuyan __ 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] lme help configuring random effects
Hi Everyone, Sorry to ask what I think is a basic question but I really haven't found my answer yet in the archives. I am trying to run a mixed effects model in R using the lme package. My experiment is such that I am interested in the effects of Temperature (2 levels) and Species (3 levels) on Growth. I collected individuals from three populations within each species. Because individuals within a population are potentially more similar to each other than individuals among populations, I want to include population as a random factor in my model. I would have thought that I would structure the model as follows: z-lme(Growth~Temp*Species, random=~1|Species/Population) But the summary for this model includes NAs (e.g. for two of the species). I've considered a model such as z--lme(Growth~Temp*Species,random=~1|Population) but it strikes me that I need to associate/nest population and species (as not all the Species treatments are applicable to every population). I'm also confused as to the naming of populations. Currently I've got the populations named 1,2,3,4,5 ...9. Should I be naming them A1,A2,A3, B1,B2,B3, C1,C2, C3 (where the letters represent different species)? When does that matter? Any help with this would be greatly appreciated! [[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] install.packages on windows
You can look for different versions of that package and try manually installing the lower version. On Fri, Sep 28, 2012 at 11:25 PM, Uwe Ligges lig...@statistik.tu-dortmund.de wrote: On 28.09.2012 00:32, Duncan Murdoch wrote: On 12-09-27 2:53 PM, Anju R wrote: Sometimes when I try to install certain packages I get a warning message. For example, I tried to install the package Imtest on windows R version 2.15.1 and got the following message: Warning message: package Imtest is not available (for R version 2.15.1) How can I install the above package? Why do I get the above Warning message? It probably means exactly what it says, except that the information is about the mirror you are using. I would try another mirror. If that doesn't solve it, then it probably means that the package is really not available for 2.15.1. You can look on the cran.r-project.org website for information about it, and probably download the source from there, but you will probably need to fix whatever is wrong with it before it will work. Or in other words: There is no such package Imtest on CRAN, perhaps you are looking for lmtest? Uwe Ligges Duncan Murdoch __** R-help@r-project.org mailing list https://stat.ethz.ch/mailman/**listinfo/r-helphttps://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/**posting-guide.htmlhttp://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-helphttps://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/** posting-guide.html 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] nlme: spatial autocorrelation on a sphere
On Sep 30, 2012, at 6:48 PM, Dan Bebber wrote: I have spatial data on a sphere (the Earth) for which I would like to run an gls model assuming that the errors are autcorrelated, i.e. including a corSpatial correlation in the model specification. In this case the distance metric should be calculated on the sphere, therefore metric = euclidean in (for example) corSpher would be incorrect. I would be grateful for help on how to write a new distance metric for the corSpatial function. I believe there are several ways that distances on a sphere can be calculated in R, for example the distMeeus function in the geosphere library. However, I have no idea how to write this into a corSpatial function. The aim is to end up with a metric = sphere option that calculates great circle distances between points using latitude and longitude. LMCTVTFY: http://cran.r-project.org/web/views/Spatial.html -- David Winsemius, MD Alameda, CA, USA __ 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] merge.zoo returns unmatched dates
On Mon, 1 Oct 2012, Vindoggy ! wrote: Sorry for the lack of reproducible data, but this seems to be a problem inherent to my dataset and I can't figure out where the issue is. I have several data frames set up as a time series with identical POSIXct date formats. If I keep the original data in data frame format and merge them using base merge- everything is perfect and everyone is happy. If I transform the data frames to zoo objects, and then do a merge.zoo- the data seem to become uncoupled from the original data. Even more unusual is that some dates in the new merged data set are prior to the original data set. I've attempted bellow to show what this looks like, and I hope someone has a suggestion as to what may be causing the problem. Here is one set of data in data.frame format head(Vup) Date Velocity_m/s 1 2010-01-21 07:42:00 1.217943 2 2010-01-21 07:43:00 1.624395 3 2010-01-21 07:44:00 1.526379 4 2010-01-21 07:45:00 1.456831 5 2010-01-21 07:46:00 1.245390 6 2010-01-21 07:47:00 1.374330 str(Vup) 'data.frame':7168 obs. of 2 variables: $ Date: POSIXct, format: 2010-01-21 07:42:00 2010-01-21 07:43:00 ... $ Velocity_m/s: num 1.22 1.62 1.53 1.46 1.25 ... And here is a second in data.frame format: head(PAS) Date PAS 1 2010-01-21 05:01:00 0.0013938 2 2010-01-21 05:02:00 0.0015331 3 2010-01-21 05:03:00 0.0016725 4 2010-01-21 05:04:00 0.0016725 5 2010-01-21 05:05:00 0.0012265 6 2010-01-21 05:06:00 0.0015889 str(PAS) 'data.frame':5520 obs. of 2 variables: $ Date : POSIXct, format: 2010-01-21 05:01:00 2010-01-21 05:02:00 ... $ PAS: num 0.00139 0.00153 0.00167 0.00167 0.00123 ... Using zoo: PASmin-zoo(as.matrix(PAS[,2]),as.POSIXct(PAS[,1],format=%Y-%m-%d %H:%M:%S,tz=UTC)) str(PASmin) ?zoo? series from 2010-01-21 05:01:00 to 2010-01-27 13:01:00 Data: num [1:5520, 1] 0.00139 0.00153 0.00167 0.00167 0.00123 ... - attr(*, dimnames)=List of 2 ..$ : NULL ..$ : chr PAS Index: POSIXct[1:5520], format: 2010-01-21 05:01:00 2010-01-21 05:02:00 2010-01-21 05:03:00 ... ADP_UPmin-zoo(as.matrix(Vup[,2]),as.POSIXct(Vup[,1], format=%Y-%m-%d %H:%M,tz=UTC)) str(ADP_UPmin) ?zoo? series from 2010-01-21 07:42:00 to 2010-01-26 20:12:00 Data: num [1:7168, 1] 1.22 1.62 1.53 1.46 1.25 ... - attr(*, dimnames)=List of 2 ..$ : NULL ..$ : chr UP_Velocity_m/s Index: POSIXct[1:7168], format: 2010-01-21 07:42:00 2010-01-21 07:43:00 2010-01-21 07:44:00 ... And if I merge the two zoo objects I get this: M-merge(ADP_UPmin,PASmin) head(M) UP_Velocity_m/s PAS 2010-01-20 21:01:00 NA 0.0013938 2010-01-20 21:02:00 NA 0.0015331 2010-01-20 21:03:00 NA 0.0016725 2010-01-20 21:04:00 NA 0.0016725 2010-01-20 21:05:00 NA 0.0012265 2010-01-20 21:06:00 NA 0.0015889 ?zoo? series from 2010-01-20 21:01:00 to 2010-01-27 05:01:00 Data: num [1:8499, 1:2] NA NA NA NA NA NA NA NA NA NA ... - attr(*, dimnames)=List of 2 ..$ : NULL ..$ : chr [1:2] UP_Velocity_m/s PAR Index: POSIXct[1:8499], format: 2010-01-20 21:01:00 2010-01-20 21:02:00 2010-01-20 21:03:00 ... For some reason I can not figure out, even though both the PAS data frame and PAS zoo object starts at 2010-01-21 05:01:00, once merged the PAS data starts a day earlier at 2010-01-20 21:01:00. The actual numeric data looks good, but both variables have no come uncoupled from the time series dates (The Velocity data is similarity uncoupled). And as stated before, doing an non-zoo merge on the data.frame data works fine. Anyone got any ideas what's going on? My guess is that you create both zoo series with time zone UTC but that the TZ attribute gets lost upon the merge. Then, the time is displayed in your systems time zone (which you haven't told us) which apparently is a couple of hours before UTC. On my system (which is in CET) I can create a series with UTC times R x - zoo(1:2, as.POSIXct(c(2012-01-01 00:00:00, +2012-01-01 01:00:00), format = %Y-%m-%d %H:%M:%S, tz = UTC)) R x 2012-01-01 00:00:00 2012-01-01 01:00:00 1 2 The times are in UTC as requested, but applying the c() method, they get dropped. See ?c.POSIXct. R time(x) [1] 2012-01-01 00:00:00 UTC 2012-01-01 01:00:00 UTC R c(time(x)) [1] 2012-01-01 01:00:00 CET 2012-01-01 02:00:00 CET Hence: R merge(x, x) x x 2012-01-01 01:00:00 1 1 2012-01-01 02:00:00 2 2 But you can set the system time in your R session to UTC which gives the desired result: R Sys.setenv(TZ = UTC) R merge(x, x) x x 2012-01-01 00:00:00 1 1 2012-01-01 01:00:00 2 2 hth, Z [[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
Re: [R] nlme: spatial autocorrelation on a sphere
Thanks, but the problem is quite specific and not addressed on the Spatial Data taskview page. Quite specifically, I would like to know how to edit corSpatial functions to calculate great circle distances. The Bayesian equivalent, georamps in the ramps package, is able to do this, therefore I imagine it must be possible for nlme. Dan From: David Winsemius [dwinsem...@comcast.net] Sent: 01 October 2012 08:38 To: Dan Bebber Cc: r-help@r-project.org Subject: Re: [R] nlme: spatial autocorrelation on a sphere On Sep 30, 2012, at 6:48 PM, Dan Bebber wrote: I have spatial data on a sphere (the Earth) for which I would like to run an gls model assuming that the errors are autcorrelated, i.e. including a corSpatial correlation in the model specification. In this case the distance metric should be calculated on the sphere, therefore metric = euclidean in (for example) corSpher would be incorrect. I would be grateful for help on how to write a new distance metric for the corSpatial function. I believe there are several ways that distances on a sphere can be calculated in R, for example the distMeeus function in the geosphere library. However, I have no idea how to write this into a corSpatial function. The aim is to end up with a metric = sphere option that calculates great circle distances between points using latitude and longitude. LMCTVTFY: http://cran.r-project.org/web/views/Spatial.html -- David Winsemius, MD Alameda, CA, USA __ 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] adjust font in ggplot2 to LaTeX document
Hi, how can i adjust the font in a ggplot2 qplot so that it will look similar to the LaTeX font? Computer Modern Sans Serif in the same size would be nice. My output device is ggsave(filename=test.pdf, width=5.5, height=3, dpi=300) and i will include the graphic with 5.5 inch in LaTeX. I found some pages about that topic but all solutions that i found have been very complicate and confusing to me. But the articles have been from 2009 too, so i guess there will be some easy solution today... Is theme the keyword to find the solution? Or will i have to include a binary font file? Can someone give me a link to an example? Kind regards and thanks a lot, -- Jonas Stein n...@jonasstein.de __ 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] merge.zoo returns unmatched dates
You nailed it, It was a tz issue. Thanks! Date: Mon, 1 Oct 2012 09:41:35 +0200 From: achim.zeil...@uibk.ac.at To: vindo...@hotmail.com CC: r-help@r-project.org Subject: Re: [R] merge.zoo returns unmatched dates On Mon, 1 Oct 2012, Vindoggy ! wrote: Sorry for the lack of reproducible data, but this seems to be a problem inherent to my dataset and I can't figure out where the issue is. I have several data frames set up as a time series with identical POSIXct date formats. If I keep the original data in data frame format and merge them using base merge- everything is perfect and everyone is happy. If I transform the data frames to zoo objects, and then do a merge.zoo- the data seem to become uncoupled from the original data. Even more unusual is that some dates in the new merged data set are prior to the original data set. I've attempted bellow to show what this looks like, and I hope someone has a suggestion as to what may be causing the problem. Here is one set of data in data.frame format head(Vup) Date Velocity_m/s 1 2010-01-21 07:42:00 1.217943 2 2010-01-21 07:43:00 1.624395 3 2010-01-21 07:44:00 1.526379 4 2010-01-21 07:45:00 1.456831 5 2010-01-21 07:46:00 1.245390 6 2010-01-21 07:47:00 1.374330 str(Vup) 'data.frame':7168 obs. of 2 variables: $ Date: POSIXct, format: 2010-01-21 07:42:00 2010-01-21 07:43:00 ... $ Velocity_m/s: num 1.22 1.62 1.53 1.46 1.25 ... And here is a second in data.frame format: head(PAS) Date PAS 1 2010-01-21 05:01:00 0.0013938 2 2010-01-21 05:02:00 0.0015331 3 2010-01-21 05:03:00 0.0016725 4 2010-01-21 05:04:00 0.0016725 5 2010-01-21 05:05:00 0.0012265 6 2010-01-21 05:06:00 0.0015889 str(PAS) 'data.frame':5520 obs. of 2 variables: $ Date : POSIXct, format: 2010-01-21 05:01:00 2010-01-21 05:02:00 ... $ PAS: num 0.00139 0.00153 0.00167 0.00167 0.00123 ... Using zoo: PASmin-zoo(as.matrix(PAS[,2]),as.POSIXct(PAS[,1],format=%Y-%m-%d %H:%M:%S,tz=UTC)) str(PASmin) ?zoo? series from 2010-01-21 05:01:00 to 2010-01-27 13:01:00 Data: num [1:5520, 1] 0.00139 0.00153 0.00167 0.00167 0.00123 ... - attr(*, dimnames)=List of 2 ..$ : NULL ..$ : chr PAS Index: POSIXct[1:5520], format: 2010-01-21 05:01:00 2010-01-21 05:02:00 2010-01-21 05:03:00 ... ADP_UPmin-zoo(as.matrix(Vup[,2]),as.POSIXct(Vup[,1], format=%Y-%m-%d %H:%M,tz=UTC)) str(ADP_UPmin) ?zoo? series from 2010-01-21 07:42:00 to 2010-01-26 20:12:00 Data: num [1:7168, 1] 1.22 1.62 1.53 1.46 1.25 ... - attr(*, dimnames)=List of 2 ..$ : NULL ..$ : chr UP_Velocity_m/s Index: POSIXct[1:7168], format: 2010-01-21 07:42:00 2010-01-21 07:43:00 2010-01-21 07:44:00 ... And if I merge the two zoo objects I get this: M-merge(ADP_UPmin,PASmin) head(M) UP_Velocity_m/s PAS 2010-01-20 21:01:00 NA 0.0013938 2010-01-20 21:02:00 NA 0.0015331 2010-01-20 21:03:00 NA 0.0016725 2010-01-20 21:04:00 NA 0.0016725 2010-01-20 21:05:00 NA 0.0012265 2010-01-20 21:06:00 NA 0.0015889 ?zoo? series from 2010-01-20 21:01:00 to 2010-01-27 05:01:00 Data: num [1:8499, 1:2] NA NA NA NA NA NA NA NA NA NA ... - attr(*, dimnames)=List of 2 ..$ : NULL ..$ : chr [1:2] UP_Velocity_m/s PAR Index: POSIXct[1:8499], format: 2010-01-20 21:01:00 2010-01-20 21:02:00 2010-01-20 21:03:00 ... For some reason I can not figure out, even though both the PAS data frame and PAS zoo object starts at 2010-01-21 05:01:00, once merged the PAS data starts a day earlier at 2010-01-20 21:01:00. The actual numeric data looks good, but both variables have no come uncoupled from the time series dates (The Velocity data is similarity uncoupled). And as stated before, doing an non-zoo merge on the data.frame data works fine. Anyone got any ideas what's going on? My guess is that you create both zoo series with time zone UTC but that the TZ attribute gets lost upon the merge. Then, the time is displayed in your systems time zone (which you haven't told us) which apparently is a couple of hours before UTC. On my system (which is in CET) I can create a series with UTC times R x - zoo(1:2, as.POSIXct(c(2012-01-01 00:00:00, +2012-01-01 01:00:00), format = %Y-%m-%d %H:%M:%S, tz = UTC)) R x 2012-01-01 00:00:00 2012-01-01 01:00:00 1 2 The times are in UTC as requested, but applying the c() method, they get dropped. See ?c.POSIXct. R time(x) [1] 2012-01-01 00:00:00 UTC 2012-01-01 01:00:00 UTC R c(time(x)) [1] 2012-01-01 01:00:00 CET 2012-01-01 02:00:00 CET Hence: R merge(x, x) x x 2012-01-01 01:00:00 1 1 2012-01-01 02:00:00 2 2 But you can set
Re: [R] REML - quasipoisson
hi the explanation seems to be that fix.family.ls is invoked to define a saturated log likelihood for quasi families, which means that the F2q term which I had thought to be undefined, is actually defined in this example as F2q--500 * log(phiq)/2 The formula then becomes F1q-F2q+F3q-F4q which does coincide with m3$gcv Note that when phiq=1, the F2q term vanishes. greg hi I'm puzzled as to the relation between the REML score computed by gam and the formula (4) on p.4 here: http://opus.bath.ac.uk/22707/1/Wood_JRSSB_2011_73_1_3.pdf I'm ok with this for poisson, or for quasipoisson when phi=1. However, when phi differs from 1, I'm stuck. #simulate some data library(mgcv) set.seed(1) x1-runif(500) x2-rnorm(500) linp--0.5+x1+exp(-x2^2/2)*sin(4*x2) y-rpois(500,exp(linp)) ##poisson #phi=1 m1-gam(y~s(x1)+s(x2),family=poisson,method=REML) phi-m1$scale #formula #1st term S1-m1$smooth[[1]]$S[[1]]*m1$sp[1] S2-m1$smooth[[2]]$S[[1]]*m1$sp[2] S-matrix(0,19,19) for (i in 2:10) { for (j in 2:10) { S[i,j]=S1[i-1,j-1] S[i+9,j+9]=S2[i-1,j-1] } } beta-m1$coef #penalised deviance Dp-m1$dev+t(beta)%*%S%*%beta F1-Dp/(2*phi) #2nd term F2-sum(ifelse(y==0,0,y*log(y)-y-log(factorial(y #3rd term X-predict(m1,type=lpmatrix) W-diag(fitted(m1)) H-t(X)%*%W%*%X ldhs-determinant(H+S,log=TRUE)$modulus[1] eigS-eigen(S,only.values=TRUE)$val lds-sum(log(eigS[1:16])) F3-(ldhs-lds)/2 #4th term Mp=3 F4-Mp/2*log(2*pi*phi) F1-F2+F3-F4 m1$gcv #reml score = formula ##quasipoisson with scale = 1 #fitting is identical to the poisson case #F1, F3, and F4 unchanged but F2 is now undefined m2-gam(y~s(x1)+s(x2),family=quasipoisson,method=REML,scale=1) F1+F3-F4 m2$gcv #reml score = formula with F2 omitted ##quasipoisson with unknown scale m3-gam(y~s(x1)+s(x2),family=quasipoisson,method=REML,scale=-1) phiq-m3$scale #1st term S1q-m3$smooth[[1]]$S[[1]]*m3$sp[1] S2q-m3$smooth[[2]]$S[[1]]*m3$sp[2] Sq-matrix(0,19,19) for (i in 2:10) { for (j in 2:10) { Sq[i,j]=S1q[i-1,j-1] Sq[i+9,j+9]=S2q[i-1,j-1] } } betaq-m3$coef #penalised deviance Dpq-m3$dev+t(betaq)%*%Sq%*%betaq F1q-Dpq/(2*phiq) #2nd term undefined #3rd term Xq-predict(m3,type=lpmatrix) Wq-diag(fitted(m3)) Hq-t(Xq)%*%Wq%*%Xq ldhsq-determinant(Hq+Sq,log=TRUE)$modulus[1] eigSq-eigen(Sq,only.values=TRUE)$val ldsq-sum(log(eigSq[1:16])) F3q-(ldhsq-ldsq)/2 #4th term Mp=3 F4q-Mp/2*log(2*pi*phiq) F1q+F3q-F4q m3$gcv #quite different #but if phiq is replaced by the Pearson estimate of the scale P-sum((y-fitted(m3))^2/fitted(m3)) phip-P/(500-sum(m3$edf)) F1p-Dpq/(2*phip) F3p-F3q #third term independent of scale F4p-Mp/2*log(2*pi*phip) F1p+F3p-F4p m3$gcv #closer but still wrong #is there a value of phi which makes this work? f-function(t) (Dpq/(2*t) + F3q + Mp/2*log(2*pi*t) - m3$gcv)^2 optimise(f,interval=c(0.5,1.5))$min-phix Dpq/(2*phix) + F3q + Mp/2*log(2*pi*phix) m3$gcv #but what is phix - not the Pearson estimate or the scale returned by gam? thanks Greg Dropkin __ 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] [Fwd: REML - quasipoisson]
Hi Greg, For quasi families I've used extended quasi-likelihood (see Mccullagh and Nelder, Generalized Linear Models 2nd ed, section 9.6) in place of the likelihood/quasi-likelihood in the expression for the (RE)ML score. I hadn't realised that this was possible before the paper was published. best, Simon ps. sorry for slow reply, the original message slipped through my filter for mgcv related stuff hi I'm puzzled as to the relation between the REML score computed by gam and the formula (4) on p.4 here: http://opus.bath.ac.uk/22707/1/Wood_JRSSB_2011_73_1_3.pdf I'm ok with this for poisson, or for quasipoisson when phi=1. However, when phi differs from 1, I'm stuck. #simulate some data library(mgcv) set.seed(1) x1-runif(500) x2-rnorm(500) linp--0.5+x1+exp(-x2^2/2)*sin(4*x2) y-rpois(500,exp(linp)) ##poisson #phi=1 m1-gam(y~s(x1)+s(x2),family=poisson,method=REML) phi-m1$scale #formula #1st term S1-m1$smooth[[1]]$S[[1]]*m1$sp[1] S2-m1$smooth[[2]]$S[[1]]*m1$sp[2] S-matrix(0,19,19) for (i in 2:10) { for (j in 2:10) { S[i,j]=S1[i-1,j-1] S[i+9,j+9]=S2[i-1,j-1] } } beta-m1$coef #penalised deviance Dp-m1$dev+t(beta)%*%S%*%beta F1-Dp/(2*phi) #2nd term F2-sum(ifelse(y==0,0,y*log(y)-y-log(factorial(y #3rd term X-predict(m1,type=lpmatrix) W-diag(fitted(m1)) H-t(X)%*%W%*%X ldhs-determinant(H+S,log=TRUE)$modulus[1] eigS-eigen(S,only.values=TRUE)$val lds-sum(log(eigS[1:16])) F3-(ldhs-lds)/2 #4th term Mp=3 F4-Mp/2*log(2*pi*phi) F1-F2+F3-F4 m1$gcv #reml score = formula ##quasipoisson with scale = 1 #fitting is identical to the poisson case #F1, F3, and F4 unchanged but F2 is now undefined m2-gam(y~s(x1)+s(x2),family=quasipoisson,method=REML,scale=1) F1+F3-F4 m2$gcv #reml score = formula with F2 omitted ##quasipoisson with unknown scale m3-gam(y~s(x1)+s(x2),family=quasipoisson,method=REML,scale=-1) phiq-m3$scale #1st term S1q-m3$smooth[[1]]$S[[1]]*m3$sp[1] S2q-m3$smooth[[2]]$S[[1]]*m3$sp[2] Sq-matrix(0,19,19) for (i in 2:10) { for (j in 2:10) { Sq[i,j]=S1q[i-1,j-1] Sq[i+9,j+9]=S2q[i-1,j-1] } } betaq-m3$coef #penalised deviance Dpq-m3$dev+t(betaq)%*%Sq%*%betaq F1q-Dpq/(2*phiq) #2nd term undefined #3rd term Xq-predict(m3,type=lpmatrix) Wq-diag(fitted(m3)) Hq-t(Xq)%*%Wq%*%Xq ldhsq-determinant(Hq+Sq,log=TRUE)$modulus[1] eigSq-eigen(Sq,only.values=TRUE)$val ldsq-sum(log(eigSq[1:16])) F3q-(ldhsq-ldsq)/2 #4th term Mp=3 F4q-Mp/2*log(2*pi*phiq) F1q+F3q-F4q m3$gcv #quite different #but if phiq is replaced by the Pearson estimate of the scale P-sum((y-fitted(m3))^2/fitted(m3)) phip-P/(500-sum(m3$edf)) F1p-Dpq/(2*phip) F3p-F3q #third term independent of scale F4p-Mp/2*log(2*pi*phip) F1p+F3p-F4p m3$gcv #closer but still wrong #is there a value of phi which makes this work? f-function(t) (Dpq/(2*t) + F3q + Mp/2*log(2*pi*t) - m3$gcv)^2 optimise(f,interval=c(0.5,1.5))$min-phix Dpq/(2*phix) + F3q + Mp/2*log(2*pi*phix) m3$gcv #but what is phix - not the Pearson estimate or the scale returned by gam? thanks Greg Dropkin -- Simon Wood, Mathematical Science, University of Bath BA2 7AY UK +44 (0)1225 386603 http://people.bath.ac.uk/sw283 __ 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] nlme: spatial autocorrelation on a sphere
On Oct 1, 2012, at 12:59 AM, Dan Bebber wrote: Thanks, but the problem is quite specific and not addressed on the Spatial Data taskview page. Quite specifically, I would like to know how to edit corSpatial functions to calculate great circle distances. The Bayesian equivalent, georamps in the ramps package, is able to do this, therefore I imagine it must be possible for nlme. Can't you use the corStruct functions in pkg::ramps? They allow specification of the 'haversine' metric. The corR* functions inherit from class corStruct. -- David. Dan From: David Winsemius [dwinsem...@comcast.net] Sent: 01 October 2012 08:38 To: Dan Bebber Cc: r-help@r-project.org Subject: Re: [R] nlme: spatial autocorrelation on a sphere On Sep 30, 2012, at 6:48 PM, Dan Bebber wrote: I have spatial data on a sphere (the Earth) for which I would like to run an gls model assuming that the errors are autcorrelated, i.e. including a corSpatial correlation in the model specification. In this case the distance metric should be calculated on the sphere, therefore metric = euclidean in (for example) corSpher would be incorrect. I would be grateful for help on how to write a new distance metric for the corSpatial function. I believe there are several ways that distances on a sphere can be calculated in R, for example the distMeeus function in the geosphere library. However, I have no idea how to write this into a corSpatial function. The aim is to end up with a metric = sphere option that calculates great circle distances between points using latitude and longitude. LMCTVTFY: http://cran.r-project.org/web/views/Spatial.html -- David Winsemius, MD Alameda, CA, USA David Winsemius, MD Alameda, CA, USA __ 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] adjust font in ggplot2 to LaTeX document
Jonas Stein news at jonasstein.de writes: Hi, how can i adjust the font in a ggplot2 qplot so that it will look similar to the LaTeX font? Computer Modern Sans Serif in the same size would be nice. My output device is ggsave(filename=test.pdf, width=5.5, height=3, dpi=300) and i will include the graphic with 5.5 inch in LaTeX. The easiest thing to do is to use the tikzDevice device for graphics output, and in turn the easiest thing to do is to use it within knitr ... google 'tikzDevice' ... __ 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] installing a package from Archive
Dear all, I'm trying to install in R the package blme that has been removed from R cran but it is still in the Archive: http://cran.r-project.org/src/contrib/Archive/blme/ I have a Mac OSX and I have been trying with the following command from the Shell: R CMD INSTALL --build http://cran.r-project.org/src/contrib/Archive/blme/blme_0.01-4.tar.gz getting the following answer: invalid package ‘0.01-4.tar.gz’ Can somebody help me? Are there other ways? Thanks in advance, Grazia -- M. Grazia Pittau, Ph.D. Associate Professor Faculty of Statistics Sapienza, University of Rome P.le Aldo Moro, 5 00185 ROMA, ITALY gra...@stat.columbia.edu grazia.pit...@uniroma1.it Phone: +39 06 49910782 Fax: +39 06 4959241 __ 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] installing a package from Archive
On 01/10/2012 7:48 AM, gra...@stat.columbia.edu wrote: Dear all, I'm trying to install in R the package blme that has been removed from R cran but it is still in the Archive: http://cran.r-project.org/src/contrib/Archive/blme/ I have a Mac OSX and I have been trying with the following command from the Shell: R CMD INSTALL --build http://cran.r-project.org/src/contrib/Archive/blme/blme_0.01-4.tar.gz getting the following answer: invalid package ‘0.01-4.tar.gz’ Can somebody help me? Are there other ways? You need a file, not a URL. So first download the file, then install it. (And you don't need the --build command line option.) Duncan Murdoch __ 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] installing a package from Archive
On 01/10/2012 12:48, gra...@stat.columbia.edu wrote: Dear all, I'm trying to install in R the package blme that has been removed from R cran but it is still in the Archive: http://cran.r-project.org/src/contrib/Archive/blme/ I have a Mac OSX and I have been trying with the following command from the Shell: R CMD INSTALL --build http://cran.r-project.org/src/contrib/Archive/blme/blme_0.01-4.tar.gz getting the following answer: invalid package ‘0.01-4.tar.gz’ Can somebody help me? Are there other ways? You have to download the file first. 'INSTALL' does not know about URLs, as reading ?INSTALL would have reminded you. Both ‘lib’ and the elements of ‘pkgs’ may be absolute or relative path names of directories. ‘pkgs’ may also contain names of package archive files: these are then extracted to a temporary directory. These are tarballs containing a single directory, optionally compressed by ‘gzip’, ‘bzip2’, ‘xz’ or ‘compress’. Finally, binary package archive files (as created by ‘R CMD INSTALL --build’) can be supplied. Thanks in advance, Grazia __ 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. -- Brian D. Ripley, rip...@stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ R-help@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] merge.zoo returns unmatched dates
HI, You can also try this: Vup-read.table(text= Date, Velocity_m/s 2010-01-21 07:42:00, 1.217943 2010-01-21 07:43:00, 1.624395 2010-01-21 07:44:00, 1.526379 2010-01-21 07:45:00, 1.456831 2010-01-21 07:46:00, 1.245390 2010-01-21 07:47:00, 1.374330 ,sep=,,header=TRUE,stringsAsFactors=FALSE) PAS-read.table(text= Date, PAS 2010-01-21 05:01:00, 0.0013938 2010-01-21 05:02:00, 0.0015331 2010-01-21 05:03:00, 0.0016725 2010-01-21 05:04:00, 0.0016725 2010-01-21 05:05:00, 0.0012265 2010-01-21 05:06:00, 0.0015889 ,sep=,,header=TRUE,stringsAsFactors=FALSE) library(xts) PAS$Date-as.POSIXct(PAS$Date,format=%Y-%m-%d %H:%M:%S,tz=UTC) Vup$Date-as.POSIXct(Vup$Date,format=%Y-%m-%d %H:%M:%S,tz=UTC) Vupxt-xts(Vup[,2],order.by=Vup[,1],tzone=UTC) PASxt-xts(PAS[,2],order.by=PAS[,1],tzone=UTC) VUPPASxt- merge(Vupxt,PASxt) VUPPASzoo-zoo(VUPPASxt) VUPPASzoo # Vupxt PASxt #2010-01-21 05:01:00 NA 0.0013938 #2010-01-21 05:02:00 NA 0.0015331 #2010-01-21 05:03:00 NA 0.0016725 #2010-01-21 05:04:00 NA 0.0016725 #2010-01-21 05:05:00 NA 0.0012265 #2010-01-21 05:06:00 NA 0.0015889 #2010-01-21 07:42:00 1.217943 NA #2010-01-21 07:43:00 1.624395 NA #2010-01-21 07:44:00 1.526379 NA #2010-01-21 07:45:00 1.456831 NA #2010-01-21 07:46:00 1.245390 NA #2010-01-21 07:47:00 1.374330 NA str(VUPPASzoo) #‘zoo’ series from 2010-01-21 05:01:00 to 2010-01-21 07:47:00 # Data: num [1:12, 1:2] NA NA NA NA NA ... #- attr(*, dimnames)=List of 2 #..$ : chr [1:12] 2010-01-21 05:01:00 2010-01-21 05:02:00 2010-01-21 05:03:00 2010-01-21 05:04:00 ... #..$ : chr [1:2] Vupxt PASxt #Index: POSIXct[1:12], format: 2010-01-21 05:01:00 2010-01-21 05:02:00 ... A.K. - Original Message - From: Vindoggy ! vindo...@hotmail.com To: r-help@r-project.org Cc: Sent: Monday, October 1, 2012 2:29 AM Subject: [R] merge.zoo returns unmatched dates Sorry for the lack of reproducible data, but this seems to be a problem inherent to my dataset and I can't figure out where the issue is. I have several data frames set up as a time series with identical POSIXct date formats. If I keep the original data in data frame format and merge them using base merge- everything is perfect and everyone is happy. If I transform the data frames to zoo objects, and then do a merge.zoo- the data seem to become uncoupled from the original data. Even more unusual is that some dates in the new merged data set are prior to the original data set. I've attempted bellow to show what this looks like, and I hope someone has a suggestion as to what may be causing the problem. Here is one set of data in data.frame format head(Vup) Date Velocity_m/s 1 2010-01-21 07:42:00 1.217943 2 2010-01-21 07:43:00 1.624395 3 2010-01-21 07:44:00 1.526379 4 2010-01-21 07:45:00 1.456831 5 2010-01-21 07:46:00 1.245390 6 2010-01-21 07:47:00 1.374330 str(Vup) 'data.frame': 7168 obs. of 2 variables: $ Date : POSIXct, format: 2010-01-21 07:42:00 2010-01-21 07:43:00 ... $ Velocity_m/s: num 1.22 1.62 1.53 1.46 1.25 ... And here is a second in data.frame format: head(PAS) Date PAS 1 2010-01-21 05:01:00 0.0013938 2 2010-01-21 05:02:00 0.0015331 3 2010-01-21 05:03:00 0.0016725 4 2010-01-21 05:04:00 0.0016725 5 2010-01-21 05:05:00 0.0012265 6 2010-01-21 05:06:00 0.0015889 str(PAS) 'data.frame': 5520 obs. of 2 variables: $ Date : POSIXct, format: 2010-01-21 05:01:00 2010-01-21 05:02:00 ... $ PAS: num 0.00139 0.00153 0.00167 0.00167 0.00123 ... Using zoo: PASmin-zoo(as.matrix(PAS[,2]),as.POSIXct(PAS[,1],format=%Y-%m-%d %H:%M:%S,tz=UTC)) str(PASmin) ‘zoo’ series from 2010-01-21 05:01:00 to 2010-01-27 13:01:00 Data: num [1:5520, 1] 0.00139 0.00153 0.00167 0.00167 0.00123 ... - attr(*, dimnames)=List of 2 ..$ : NULL ..$ : chr PAS Index: POSIXct[1:5520], format: 2010-01-21 05:01:00 2010-01-21 05:02:00 2010-01-21 05:03:00 ... ADP_UPmin-zoo(as.matrix(Vup[,2]),as.POSIXct(Vup[,1], format=%Y-%m-%d %H:%M,tz=UTC)) str(ADP_UPmin) ‘zoo’ series from 2010-01-21 07:42:00 to 2010-01-26 20:12:00 Data: num [1:7168, 1] 1.22 1.62 1.53 1.46 1.25 ... - attr(*, dimnames)=List of 2 ..$ : NULL ..$ : chr UP_Velocity_m/s Index: POSIXct[1:7168], format: 2010-01-21 07:42:00 2010-01-21 07:43:00 2010-01-21 07:44:00 ... And if I merge the two zoo objects I get this: M-merge(ADP_UPmin,PASmin) head(M) UP_Velocity_m/s PAS 2010-01-20 21:01:00 NA 0.0013938 2010-01-20 21:02:00 NA 0.0015331 2010-01-20 21:03:00 NA 0.0016725 2010-01-20 21:04:00 NA 0.0016725 2010-01-20 21:05:00 NA 0.0012265 2010-01-20 21:06:00 NA 0.0015889 ‘zoo’ series from 2010-01-20 21:01:00 to 2010-01-27 05:01:00 Data: num [1:8499, 1:2] NA NA NA NA
Re: [R] lme help configuring random effects
Julie Lee-Yaw julleeyaw at yahoo.ca writes: [snip] I am trying to run a mixed effects model in R using the lme package. My experiment is such that I am interested in the effects of Temperature (2 levels) and Species (3 levels) on Growth. I collected individuals from three populations within each species. Because individuals within a population are potentially more similar to each other than individuals among populations, I want to include population as a random factor in my model. You may have some practical difficulties incorporating a random effect with only three populations ... I would have thought that I would structure the model as follows: z-lme(Growth~Temp*Species, random=~1|Species/Population) But the summary for this model includes NAs (e.g. for two of the species). I've considered a model such as z- lme(Growth~Temp*Species,random=~1|Population) I think you need z - lme(Growth~Temp*Species, random=~1|Species:Population) Your specification (Species/Population) expands to Species+ Species:Population , which ends up including Species as both a random and a fixed factor. You probably don't have the data, but in principle you should consider random=~Temp|Species:Population (see a paper by Schielzeth on incorporating treatment-by-block interactions) [snip] I'm also confused as to the naming of populations. Currently I've got the populations named 1,2,3,4,5 ...9. Should I be naming them A1,A2,A3, B1,B2,B3, C1,C2, C3 (where the letters represent different species)? When does that matter? This is the difference between implicit and explicit nesting. In general naming them A1, ... C1, ... is better, because it reduces the probability of making mistakes. http://glmm.wikidot.com/faq may be useful. Further questions should probably go to the r-sig-mixed-models mailing list. __ 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] svyboxplot - library (survey)
Hello, I have used the library (survey) package for boxplots using the following code. Could anyone please tell me why I am getting only 1 boxplot instead of 2 boxplots (1-SPD, 2-No SPD). What changes in the following code would be required to get 2 boxplots in the same plot frame? Thanks, Pradip ### nhis - svydesign (id=~psu, strat=~stratum, weights=~wt8, data=tor, nest=TRUE) svyboxplot (dthage~xspd2, subset (nhis, mortstat==1), col=gray80, varwidth=TRUE, ylab=Age at Death, xlab=SPD Status: 1-SPD, 2=No SPD) Pradip K. Muhuri Statistician Substance Abuse Mental Health Services Administration The Center for Behavioral Health Statistics and Quality Division of Population Surveys 1 Choke Cherry Road, Room 2-1071 Rockville, MD 20857 Tel: 240-276-1070 Fax: 240-276-1260 e-mail: pradip.muh...@samhsa.hhs.gov The Center for Behavioral Health Statistics and Quality your feedback. Please click on the following link to complete a brief customer survey: http://cbhsqsurvey.samhsa.gov vide 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] installing a package from Archive
On Mon, 1 Oct 2012, Duncan Murdoch wrote: Thanks a lot. Following your suggestions, I have done the following: zercona:Downloads graziapittau$ R CMD INSTALL blme 0.01-4.tar.gz Warning: invalid package ‘0.01-4.tar.gz’ * installing to library ‘/Library/Frameworks/R.framework/Versions/2.15/Resources/library’ * installing *source* package ‘blme’ ... ** package ‘blme’ successfully unpacked and MD5 sums checked ** libs *** arch - i386 sh: make: command not found ERROR: compilation failed for package ‘blme’ * removing ‘/Library/Frameworks/R.framework/Versions/2.15/Resources/library/blme’ * restoring previous ‘/Library/Frameworks/R.framework/Versions/2.15/Resources/library/blme’ Unfortunately, I still get an error message and from Rthe following: library(blme) Error in library(blme) : ‘blme’ is not a valid installed package Where am I doing wrong? Thanks a lot. Grazia On 01/10/2012 7:48 AM, gra...@stat.columbia.edu wrote: Dear all, I'm trying to install in R the package blme that has been removed from R cran but it is still in the Archive: http://cran.r-project.org/src/contrib/Archive/blme/ I have a Mac OSX and I have been trying with the following command from the Shell: R CMD INSTALL --build http://cran.r-project.org/src/contrib/Archive/blme/blme_0.01-4.tar.gz getting the following answer: invalid package ‘0.01-4.tar.gz’ Can somebody help me? Are there other ways? You need a file, not a URL. So first download the file, then install it. (And you don't need the --build command line option.) Duncan Murdoch -- M. Grazia Pittau, Ph.D. Associate Professor Faculty of Statistics Sapienza, University of Rome P.le Aldo Moro, 5 00185 ROMA, ITALY gra...@stat.columbia.edu grazia.pit...@uniroma1.it Phone: +39 06 49910782 Fax: +39 06 4959241 __ 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] installing a package from Archive
On 01/10/2012 10:02 AM, gra...@stat.columbia.edu wrote: On Mon, 1 Oct 2012, Duncan Murdoch wrote: Thanks a lot. Following your suggestions, I have done the following: zercona:Downloads graziapittau$ R CMD INSTALL blme 0.01-4.tar.gz Warning: invalid package ‘0.01-4.tar.gz’ The filename should not have a space in it. I think the original had an underscore; if that's been lost, you need to rename the file to something with no spaces. (You can probably quote the filename, but the details of how to do that depend on your shell; it's easier just to avoid filenames that contain spaces.) Duncan Murdoch * installing to library ‘/Library/Frameworks/R.framework/Versions/2.15/Resources/library’ * installing *source* package ‘blme’ ... ** package ‘blme’ successfully unpacked and MD5 sums checked ** libs *** arch - i386 sh: make: command not found ERROR: compilation failed for package ‘blme’ * removing ‘/Library/Frameworks/R.framework/Versions/2.15/Resources/library/blme’ * restoring previous ‘/Library/Frameworks/R.framework/Versions/2.15/Resources/library/blme’ Unfortunately, I still get an error message and from Rthe following: library(blme) Error in library(blme) : ‘blme’ is not a valid installed package Where am I doing wrong? Thanks a lot. Grazia On 01/10/2012 7:48 AM, gra...@stat.columbia.edu wrote: Dear all, I'm trying to install in R the package blme that has been removed from R cran but it is still in the Archive: http://cran.r-project.org/src/contrib/Archive/blme/ I have a Mac OSX and I have been trying with the following command from the Shell: R CMD INSTALL --build http://cran.r-project.org/src/contrib/Archive/blme/blme_0.01-4.tar.gz getting the following answer: invalid package ‘0.01-4.tar.gz’ Can somebody help me? Are there other ways? You need a file, not a URL. So first download the file, then install it. (And you don't need the --build command line option.) Duncan Murdoch __ 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] svyboxplot - library (survey)
using a slight modification of the example shown in ?svyboxplot # load survey library library(survey) # load example data data(api) # create an example svydesign dstrat - svydesign(id = ~1, strata = ~stype, weights = ~pw, data = apistrat, fpc = ~fpc) # set the plot window to display 1 plot x 2 plots par(mfrow=c(1,2)) # generate two example boxplots svyboxplot(enroll~stype,dstrat,all.outliers=TRUE) svyboxplot(enroll~1,dstrat) # done # alternative: not as nice # set the plot window to display 2 plots x 1 plot par(mfrow=c(2,1)) # generate two example boxplots svyboxplot(enroll~stype,dstrat,all.outliers=TRUE) svyboxplot(enroll~1,dstrat) # done On Mon, Oct 1, 2012 at 9:50 AM, Muhuri, Pradip (SAMHSA/CBHSQ) pradip.muh...@samhsa.hhs.gov wrote: Hello, I have used the library (survey) package for boxplots using the following code. Could anyone please tell me why I am getting only 1 boxplot instead of 2 boxplots (1-SPD, 2-No SPD). What changes in the following code would be required to get 2 boxplots in the same plot frame? Thanks, Pradip ### nhis - svydesign (id=~psu, strat=~stratum, weights=~wt8, data=tor, nest=TRUE) svyboxplot (dthage~xspd2, subset (nhis, mortstat==1), col=gray80, varwidth=TRUE, ylab=Age at Death, xlab=SPD Status: 1-SPD, 2=No SPD) Pradip K. Muhuri Statistician Substance Abuse Mental Health Services Administration The Center for Behavioral Health Statistics and Quality Division of Population Surveys 1 Choke Cherry Road, Room 2-1071 Rockville, MD 20857 Tel: 240-276-1070 Fax: 240-276-1260 e-mail: pradip.muh...@samhsa.hhs.gov The Center for Behavioral Health Statistics and Quality your feedback. Please click on the following link to complete a brief customer survey: http://cbhsqsurvey.samhsa.gov vide 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. [[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] adjust font in ggplot2 to LaTeX document
On 9/30/2012 11:38 AM, Jonas Stein wrote: Hi, how can i adjust the font in a ggplot2 qplot so that it will look similar to the LaTeX font? Computer Modern Sans Serif in the same size would be nice. My output device is ggsave(filename=test.pdf, width=5.5, height=3, dpi=300) and i will include the graphic with 5.5 inch in LaTeX. I found some pages about that topic but all solutions that i found have been very complicate and confusing to me. But the articles have been from 2009 too, so i guess there will be some easy solution today... There is a recent (late last month) blog post on this topic: http://blog.revolutionanalytics.com/2012/09/how-to-use-your-favorite-fonts-in-r-charts.html Is theme the keyword to find the solution? In ggplot2, the font family is controlled by the theme. font family will probably get more specific results. Or will i have to include a binary font file? Possibly, but probably not since you would be using a font already present in the (final) PDF. Can someone give me a link to an example? Kind regards and thanks a lot, -- Brian S. Diggs, PhD Senior Research Associate, Department of Surgery Oregon Health Science University __ 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] Better way of Grouping?
Thank you Jeff, The main purpose I am looking to use these subsets for is to do comparisons between groups and/or timepoints within groups. For example the difference in means between 3 and 4, or the percent difference between group L an D within group 3. The code I have provided is almost exactly as used although David accurately noted I mistakenly typed '3' when I intended to put '4'. I am just looking for some general guidance in improving how I go about my code. Because there are so many subgroups that can be formed with three layers of groups, it simply becomes a clutter and I wanted to learn if there was a 'better' way to organize groups for comparisons. Does that help clarify? Thanks again, Charles On Fri, Sep 28, 2012 at 4:25 PM, Jeff Newmiller jdnew...@dcn.davis.ca.uswrote: You have not specified the objective function you are trying to optimize with your term efficient, or what you do with all of these subsets once you have them. For notational simplification and completeness of coverage (not necessarily computational speedup) you might want to look at tapply or ddply/dlply from the plyr package. If you build lists of subsets you can index into them according to grouping value. You can use expand.grid to build all permutations of grouping values to use as indexes into those lists of subsets. To reiterate, you have not indicated what you want to do with these subsets, so there could be special-purpose functions that do what you want. As always, reproducible code leads to reproducible answers. :) --- 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. Charles Determan Jr deter...@umn.edu wrote: Hello R users, This is more of a convenience question that I hope others might find useful if there is a better answer. I work with large datasets that requires multiple parsing stages for different analysis. For example, compare group 3 vs. group 4. A more complicated comparison would be time B in group 3 of group L with B in group 4 of group L. I normally subset each group with the following type of code. data=read(...) #L v D L=data[LvD %in% c(L),] D=data[LvD %in% c(D),] #Groups 3 and 4 within L and D group3L=L[group %in% c(3),] group4L=L[group %in% c(3),] group3D=D[group %in% c(3),] group4D=D[group %in% c(3),] #Times B, S45, FR2, FR8 you get the idea Is there a more efficient way to subset groups? Thanks for any insight. Regards, Charles [[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] which() with multiple conditions
I hope someone can point me in the right direction please. I have a data frame with a column containing names. I want to identify the columns that contain names in a list. namestofind - c('fred','bill',a long list) If I only wanted to identify a single name I would use which(z$name == 'bill') What syntax would I use to identify all the rows that contain any of the names in namestofind? Thanks in advance for the pointer -- View this message in context: http://r.789695.n4.nabble.com/which-with-multiple-conditions-tp4644677.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] can't create png graphics under linux anymore
Hi, since my update of Linux from opensuse 11.4 to 12.1 and the update of R to 2.15.1 I can't produce graphics with png() or similar (like tiff, jpeg etc.) anymore. The following error occurs: png() Error in X11(paste(png::, filename, sep = ), g$width, g$height, pointsize, : unable to start Device PNG Additionally: Warning: In png() : png isn't supported of this R version (error message translated from german) Before the updates everything worked fine and still works under Windows 7 for R 2.15.1 on the same machine. I couldn't find something useful with google. As far as I know R still supports png?! Maybe the problem is related to suse? Thanks for any help, Tobias __ 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] enquiry
hi, I am new to using R ,I have 2 datasets with dates common in them ,how can i take out the common dates within them. karan [[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] calculating probability from the density function
Hello, I have a data x with normal (or very close to normal) distribution, I can plot a density distribution with density(x,...). My question is is there any way to calculate an area under this distribution (=probability) for particular range of x values, let's say for x from 0 to 2? I was not able to find any kind of simple procedure to do this. Thanks in advance for your help, Evgeny. [[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] which() with multiple conditions
On Monday, October 1, 2012, pdb wrote: I hope someone can point me in the right direction please. I have a data frame with a column containing names. I want to identify the columns that contain names in a list. namestofind - c('fred','bill',a long list) If I only wanted to identify a single name I would use which(z$name == 'bill') What syntax would I use to identify all the rows that contain any of the names in namestofind? %in% Michael Thanks in advance for the pointer -- View this message in context: http://r.789695.n4.nabble.com/which-with-multiple-conditions-tp4644677.html Sent from the R help mailing list archive at Nabble.com. __ R-help@r-project.org javascript:; 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] enquiry
On Monday, October 1, 2012, Karan Anand wrote: hi, I am new to using R ,I have 2 datasets with dates common in them ,how can i take out the common dates within them. Depending on what you mean by 'out' either merge() or setdiff() RMW karan [[alternative HTML version deleted]] __ R-help@r-project.org javascript:; 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] Retrieve hypergeometric results in large scale
I'm going to use dhyper(x, m, n, k) to get a 95% coverage. Let me use an example to explain my problem: Suppose I have a urn containing 90 red and 10 black balls. Now I wanna remove 3 from the urn. By the following codes: m-90;n-10;k-3; x-0:3 dhyper(x,m,n,k) I can obtain the probability that 0,1,2,3 red balls will be removed. 0.000742115 0.025046382 0.247680891 0.726530612 So 95% time, 2 to 3 red balls will be removed and the resultant composition will be changed to 87:10 or 88:9, the original percent of red balls will be changed from 90 to 89.69 to 90.72 then. If now I have 50:50 and again to remove 3 balls, I will obtain the probability as: 0.1212 0.3788 0.3788 0.1212 To get the resultant range of red balls for 95% time, this time all the four cases have to consider and so the resultant change of red balls will become 48.45 to 51.54 So my problem is, is there any convenient built-in function that helps extract this 95% confidence interval-like data? -- View this message in context: http://r.789695.n4.nabble.com/Retrieve-hypergeometric-results-in-large-scale-tp4644683.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] Retrieve hypergeometric results in large scale
Perhaps you should read ?dhyper and if you have a hard time parsing that, then read ?Distributions and then go back to ?dhyper --- 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. jas4710 wata...@post.com wrote: I'm going to use dhyper(x, m, n, k) to get a 95% coverage. Let me use an example to explain my problem: Suppose I have a urn containing 90 red and 10 black balls. Now I wanna remove 3 from the urn. By the following codes: m-90;n-10;k-3; x-0:3 dhyper(x,m,n,k) I can obtain the probability that 0,1,2,3 red balls will be removed. 0.000742115 0.025046382 0.247680891 0.726530612 So 95% time, 2 to 3 red balls will be removed and the resultant composition will be changed to 87:10 or 88:9, the original percent of red balls will be changed from 90 to 89.69 to 90.72 then. If now I have 50:50 and again to remove 3 balls, I will obtain the probability as: 0.1212 0.3788 0.3788 0.1212 To get the resultant range of red balls for 95% time, this time all the four cases have to consider and so the resultant change of red balls will become 48.45 to 51.54 So my problem is, is there any convenient built-in function that helps extract this 95% confidence interval-like data? -- View this message in context: http://r.789695.n4.nabble.com/Retrieve-hypergeometric-results-in-large-scale-tp4644683.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] Retrieve hypergeometric results in large scale
Homework? There's a no homework policy on this list. -- Bert On Mon, Oct 1, 2012 at 8:10 AM, Jeff Newmiller jdnew...@dcn.davis.ca.us wrote: Perhaps you should read ?dhyper and if you have a hard time parsing that, then read ?Distributions and then go back to ?dhyper --- 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. jas4710 wata...@post.com wrote: I'm going to use dhyper(x, m, n, k) to get a 95% coverage. Let me use an example to explain my problem: Suppose I have a urn containing 90 red and 10 black balls. Now I wanna remove 3 from the urn. By the following codes: m-90;n-10;k-3; x-0:3 dhyper(x,m,n,k) I can obtain the probability that 0,1,2,3 red balls will be removed. 0.000742115 0.025046382 0.247680891 0.726530612 So 95% time, 2 to 3 red balls will be removed and the resultant composition will be changed to 87:10 or 88:9, the original percent of red balls will be changed from 90 to 89.69 to 90.72 then. If now I have 50:50 and again to remove 3 balls, I will obtain the probability as: 0.1212 0.3788 0.3788 0.1212 To get the resultant range of red balls for 95% time, this time all the four cases have to consider and so the resultant change of red balls will become 48.45 to 51.54 So my problem is, is there any convenient built-in function that helps extract this 95% confidence interval-like data? -- View this message in context: http://r.789695.n4.nabble.com/Retrieve-hypergeometric-results-in-large-scale-tp4644683.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. -- 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.
[R] Problem with nls regression fit
Hi all, I got following problem in fitting the data. Any kind of suggestions are welcome beta - 3.5 d - seq(0.1,62.5,0.1) y - exp(-beta*d) y1 - y x - read.table(epidist.txt, header = TRUE) data.nls - as.data.frame(cbind(y1,x)) #attach(data.nls) nls.fit - nls(y1~dist,data.nls) Error in cll[[1L]] : object of type 'symbol' is not subsettable Best [[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] svyboxplot - library (survey)
Dear Anthony, Yes, I can follow the example code you have given. But, do you know from the code shown below (following Thomas Lumley's Complex Surveys) why I am getting the boxplot of dthage for just xspd=1, not xspd2=2? My intent is the make this code work so that I can generate similar plots on other continuous variable. Any help will be appreciated. Thanks, Pradip nhis - svydesign (id=~psu, strat=~stratum, weights=~wt8, data=tor, nest=TRUE) svyboxplot (dthage~xspd2, subset (nhis, mortstat==1), col=gray80, varwidth=TRUE, ylab=Age at Death, xlab=SPD Status: 1-SPD, 2=No SPD) Pradip K. Muhuri, PhD Statistician Substance Abuse Mental Health Services Administration The Center for Behavioral Health Statistics and Quality Division of Population Surveys 1 Choke Cherry Road, Room 2-1071 Rockville, MD 20857 Tel: 240-276-1070 Fax: 240-276-1260 e-mail: pradip.muh...@samhsa.hhs.govmailto:pradip.muh...@samhsa.hhs.gov The Center for Behavioral Health Statistics and Quality your feedback. Please click on the following link to complete a brief customer survey: http://cbhsqsurvey.samhsa.govhttp://cbhsqsurvey.samhsa.gov/ From: Anthony Damico [mailto:ajdam...@gmail.com] Sent: Monday, October 01, 2012 10:07 AM To: Muhuri, Pradip (SAMHSA/CBHSQ) Cc: R help Subject: Re: [R] svyboxplot - library (survey) using a slight modification of the example shown in ?svyboxplot # load survey library library(survey) # load example data data(api) # create an example svydesign dstrat - svydesign(id = ~1, strata = ~stype, weights = ~pw, data = apistrat, fpc = ~fpc) # set the plot window to display 1 plot x 2 plots par(mfrow=c(1,2)) # generate two example boxplots svyboxplot(enroll~stype,dstrat,all.outliers=TRUE) svyboxplot(enroll~1,dstrat) # done # alternative: not as nice # set the plot window to display 2 plots x 1 plot par(mfrow=c(2,1)) # generate two example boxplots svyboxplot(enroll~stype,dstrat,all.outliers=TRUE) svyboxplot(enroll~1,dstrat) # done On Mon, Oct 1, 2012 at 9:50 AM, Muhuri, Pradip (SAMHSA/CBHSQ) pradip.muh...@samhsa.hhs.govmailto:pradip.muh...@samhsa.hhs.gov wrote: Hello, I have used the library (survey) package for boxplots using the following code. Could anyone please tell me why I am getting only 1 boxplot instead of 2 boxplots (1-SPD, 2-No SPD). What changes in the following code would be required to get 2 boxplots in the same plot frame? Thanks, Pradip ### nhis - svydesign (id=~psu, strat=~stratum, weights=~wt8, data=tor, nest=TRUE) svyboxplot (dthage~xspd2, subset (nhis, mortstat==1), col=gray80, varwidth=TRUE, ylab=Age at Death, xlab=SPD Status: 1-SPD, 2=No SPD) Pradip K. Muhuri Statistician Substance Abuse Mental Health Services Administration The Center for Behavioral Health Statistics and Quality Division of Population Surveys 1 Choke Cherry Road, Room 2-1071 Rockville, MD 20857 Tel: 240-276-1070 Fax: 240-276-1260 e-mail: pradip.muh...@samhsa.hhs.govmailto:pradip.muh...@samhsa.hhs.gov The Center for Behavioral Health Statistics and Quality your feedback. Please click on the following link to complete a brief customer survey: http://cbhsqsurvey.samhsa.gov vide commented, minimal, self-contained, reproducible code. __ R-help@r-project.orgmailto: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] calculating probability from the density function
Forget density(). Smooth the ecdf instead. ?lowess ?smooth.spline ?monotone.smooth (in package fda, for monotone smoothing, which may be preferable) ?locfit (in package locfit) ... and many many others -- Bert On Mon, Oct 1, 2012 at 6:46 AM, Eugene Kanshin kanshin...@gmail.com wrote: Hello, I have a data x with normal (or very close to normal) distribution, I can plot a density distribution with density(x,...). My question is is there any way to calculate an area under this distribution (=probability) for particular range of x values, let's say for x from 0 to 2? I was not able to find any kind of simple procedure to do this. Thanks in advance for your help, Evgeny. [[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. -- 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] can't create png graphics under linux anymore
See the 'R Installation and Administration Manual'. And yes, it is related to the actions of whoever installed R, so please discuss it with them. On 01/10/2012 14:54, Tobias Pilz wrote: Hi, since my update of Linux from opensuse 11.4 to 12.1 and the update of R to 2.15.1 I can't produce graphics with png() or similar (like tiff, jpeg etc.) anymore. The following error occurs: png() Error in X11(paste(png::, filename, sep = ), g$width, g$height, pointsize, : unable to start Device PNG Additionally: Warning: In png() : png isn't supported of this R version (error message translated from german) Before the updates everything worked fine and still works under Windows 7 for R 2.15.1 on the same machine. I couldn't find something useful with google. As far as I know R still supports png?! Maybe the problem is related to suse? Thanks for any help, Tobias __ 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. -- Brian D. Ripley, rip...@stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] calculating probability from the density function
On Oct 1, 2012, at 6:46 AM, Eugene Kanshin wrote: Hello, I have a data x with normal (or very close to normal) distribution, I can plot a density distribution with density(x,...). My question is is there any way to calculate an area under this distribution (=probability) for particular range of x values, let's say for x from 0 to 2? I was not able to find any kind of simple procedure to do this. ? pnorm -- David Winsemius, MD Alameda, CA, USA __ 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] Hmisc describe error
Describe fails for me with a message similar to what was an issue in 2008 and got fixed according to posts. R version 2.15.0 (2012-03-30) Copyright (C) 2012 The R Foundation for Statistical Computing ISBN 3-900051-07-0 Platform: i386-pc-mingw32/i386 (32-bit) # output truncated options(chmhelp = FALSE, help_type = text) .help.ESS - help options(STERM='iESS', editor='gnuclient.exe') load('prostate.sav') library(rms) Loading required package: Hmisc Loading required package: survival Loading required package: splines Hmisc library by Frank E Harrell Jr Type library(help='Hmisc'), ?Overview, or ?Hmisc.Overview') to see overall documentation. NOTE:Hmisc no longer redefines [.factor to drop unused levels when subsetting. To get the old behavior of Hmisc type dropUnusedLevels(). Attaching package: 'Hmisc' The following object(s) are masked from 'package:survival': untangle.specials The following object(s) are masked from 'package:base': format.pval, round.POSIXt, trunc.POSIXt, units Attaching package: 'rms' The following object(s) are masked from 'package:survival': Surv describe(prostate) Error in format(dates(x)) : could not find function dates Thanks everybody. Stephen __ 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] nlme: spatial autocorrelation on a sphere
On Oct 1, 2012, at 8:34 AM, Spencer Graves wrote: On 10/1/2012 12:38 AM, David Winsemius wrote: snip LMCTVTFY: http://cran.r-project.org/web/views/Spatial.html What's LMCTVIFY? Please accept my apologies for attempting to be cute and I also apologize to Mr Bebber for reading his request far too quickly. I was trying to use (C)RAN (T)ask (V)iews as a verb. LMGTFY for LMCTVIY produced only LMGTFY, at least for me. There's space for it on Wikipedia (RTFM) after LMGTFY (which I found using LMGTFY Wikipedia). I don't think it warrants being enshrined. I might also have written: require(sos); findFn(metric spherical latitude longitude) might produce. In this instance that approach found the ramps::corRSpher function, which has a haversine metric, much more quickly than did my subsequent efforts with RSeek, which I conducted after I realized that Bebber had already made a good faith effort at identifying resources. -- David Winsemius, MD Alameda, CA, USA __ 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] Retrieve hypergeometric results in large scale
Thanks Jeff The documentation pages, if I haven't missed any crucial points, illustrate how to get probability and cumulative probability values. I can first retrieve the data structures and use Perl (I don't know how to use R...) to sort the derived ratios and sum the probability values until the cumulative probability exceeds 95%. Well, I just don't know whether such seemingly routine procedures have already been implemented... Thanks again! -- View this message in context: http://r.789695.n4.nabble.com/Retrieve-95-coverage-of-results-from-a-hypergeometric-distribution-tp4644683p4644701.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] Error messages when attempting to calculate polychoric correlation matrices
Dear R users, I am a psychology postgraduate student who is relatively new to using R. I am currently developing a psychometric scale and have run into a few problems when using R to calculate a polychoric correlation matrix for my dataset. I am trying to produce a polychoric correlation matrix for calculating ordinal reliability estimates (eg. Alpha, omega).The set consists of 439 observations of 18 variables. Each of the 18 variables was measured using a 6 point Likert scale but some variables only feature responses across 5 ordinal steps because there were no observed responses for one extreme of the item. The dataset has no missing data following multiple imputation and preliminary analyses were conducted using SPSS. I have included a command that will reproduce the dataset below, apologies for the size of the command: data.matrix18 - structure(list(SABAS02 = structure(c(6, 2, 5, 6, 4, 6, 5, 3, 5, 1, 4, 4, 6, 5, 5, 5, 4, 4, 4, 3, 4, 6, 5, 2, 5, 4, 5, 4, 4, 5, 5, 5, 2, 5, 4, 5, 5, 5, 6, 3, 6, 4, 5, 5, 5, 5, 4, 4, 5, 5, 5, 4, 5, 5, 4, 6, 6, 5, 6, 4, 4, 5, 4, 5, 5, 5, 6, 4, 4, 5, 4, 6, 5, 6, 5, 6, 6, 5, 4, 6, 6, 5, 3, 4, 2, 4, 5, 5, 6, 2, 6, 5, 4, 6, 2, 5, 3, 5, 4, 4, 5, 4, 4, 3, 4, 4, 6, 5, 6, 6, 4, 4, 5, 5, 5, 6, 5, 5, 5, 5, 4, 5, 5, 5, 4, 4, 5, 6, 4, 5, 5, 4, 5, 5, 4, 5, 5, 4, 5, 5, 5, 4, 2, 5, 4, 5, 5, 5, 5, 5, 5, 4, 5, 4, 4, 5, 4, 6, 2, 4, 6, 4, 4, 5, 4, 4, 6, 6, 5, 3, 3, 2, 4, 5, 5, 5, 5, 6, 2, 5, 5, 4, 6, 4, 5, 6, 5, 4, 5, 5, 4, 4, 5, 4, 5, 5, 3, 5, 4, 3, 5, 4, 5, 5, 5, 5, 5, 5, 5, 5, 4, 4, 5, 4, 4, 5, 5, 6, 4, 3, 5, 6, 5, 4, 5, 5, 2, 5, 5, 4, 2, 4, 5, 5, 5, 5, 5, 3, 4, 5, 2, 4, 3, 4, 5, 5, 5, 5, 5, 5, 5, 6, 6, 5, 4, 5, 6, 4, 5, 3, 5, 6, 5, 4, 5, 5, 5, 5, 5, 4, 4, 5, 6, 5, 4, 5, 5, 6, 6, 6, 4, 3, 4, 4, 5, 6, 4, 5, 5, 4, 6, 5, 5, 4, 5, 4, 4, 4, 4, 5, 3, 5, 5, 4, 3, 4, 5, 4, 4, 5, 4, 5, 5, 3, 5, 6, 5, 5, 5, 5, 5, 4, 3, 5, 4, 4, 5, 5, 6, 4, 5, 4, 2, 6, 5, 5, 5, 5, 5, 5, 5, 6, 6, 5, 5, 5, 4, 5, 5, 6, 5, 4, 4, 5, 4, 4, 5, 5, 4, 6, 4, 4, 5, 6, 5, 5, 4, 5, 4, 5, 4, 6, 5, 5, 5, 4, 5, 6, 5, 5, 5, 5, 5, 4, 3, 4, 5, 4, 4, 6, 4, 4, 4, 2, 5, 4, 4, 6, 2, 6, 4, 5, 5, 6, 5, 6, 4, 4, 3, 4, 5, 5, 6, 5, 5, 5, 5, 4, 5, 3, 5, 5, 5, 5, 6, 3, 5, 6, 5, 6, 4, 5, 4, 4, 5, 6, 2, 5, 6), value.labels = structure(c(6, 5, 4, 3, 2, 1), .Names = c(Strongly Agree, Agree, Slightly Agree, Slightly Disagree, Disagree, Strongly Disagree))), SABAS06 = structure(c(6, 6, 6, 6, 6, 6, 5, 5, 6, 3, 6, 6, 6, 5, 5, 6, 5, 6, 6, 5, 5, 6, 6, 6, 6, 5, 6, 3, 6, 5, 5, 5, 5, 6, 6, 6, 3, 6, 6, 5, 5, 4, 5, 6, 6, 5, 6, 5, 6, 6, 6, 6, 6, 6, 5, 5, 6, 6, 6, 6, 6, 6, 5, 6, 6, 6, 6, 6, 4, 6, 5, 5, 6, 6, 6, 6, 6, 6, 3, 5, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 5, 5, 6, 6, 6, 5, 5, 5, 6, 5, 6, 5, 6, 6, 6, 6, 5, 5, 6, 5, 5, 6, 2, 5, 5, 6, 5, 5, 5, 6, 5, 6, 5, 5, 6, 6, 6, 5, 5, 5, 6, 6, 5, 6, 6, 6, 4, 6, 6, 5, 4, 6, 5, 5, 6, 5, 5, 5, 6, 6, 6, 6, 5, 5, 6, 6, 5, 6, 4, 6, 6, 6, 6, 5, 5, 6, 4, 5, 6, 5, 6, 5, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 5, 5, 6, 6, 5, 5, 6, 6, 5, 6, 5, 6, 5, 4, 6, 5, 6, 6, 5, 6, 5, 6, 5, 5, 5, 5, 6, 6, 5, 6, 5, 4, 6, 5, 6, 6, 5, 5, 6, 5, 5, 6, 6, 6, 6, 5, 5, 6, 6, 5, 5, 6, 6, 4, 5, 5, 6, 6, 6, 6, 4, 6, 6, 6, 6, 6, 5, 6, 6, 6, 6, 5, 6, 6, 6, 5, 5, 6, 5, 6, 6, 6, 5, 2, 6, 6, 5, 5, 6, 5, 6, 6, 5, 6, 5, 6, 5, 6, 6, 6, 5, 5, 4, 5, 6, 6, 5, 5, 4, 4, 6, 6, 5, 4, 5, 5, 6, 6, 6, 5, 4, 5, 5, 3, 5, 5, 6, 6, 5, 5, 5, 5, 5, 5, 5, 6, 6, 6, 6, 5, 5, 5, 5, 5, 5, 4, 5, 5, 6, 5, 6, 6, 6, 5, 6, 5, 6, 6, 6, 5, 5, 6, 5, 6, 6, 6, 6, 5, 6, 6, 5, 6, 5, 6, 6, 5, 5, 5, 6, 5, 5, 5, 5, 6, 5, 6, 6, 5, 5, 6, 5, 5, 6, 6, 4, 3, 3, 4, 3, 5, 4, 5, 5, 6, 6, 6, 6, 6, 6, 2, 5, 6, 5, 6, 6, 6, 6, 6, 6, 5, 6, 5, 6, 6, 6, 6, 5, 6, 5, 5, 5, 6, 6, 6, 6, 5, 6, 5, 6, 5, 6, 6, 6, 5, 4, 6, 6, 6, 6, 6, 5), value.labels = structure(c(6, 5, 4, 3, 2, 1), .Names = c(Strongly Agree, Agree, Slightly Agree, Slightly Disagree, Disagree, Strongly Disagree))), SABAS09 = structure(c(5, 3, 5, 5, 4, 5, 4, 2, 5, 1, 4, 5, 5, 5, 3, 5, 4, 4, 6, 4, 4, 6, 5, 6, 4, 4, 4, 4, 4, 4, 5, 5, 3, 4, 5, 5, 2, 2, 5, 4, 5, 4, 6, 5, 6, 4, 3, 6, 5, 5, 5, 5, 6, 5, 4, 5, 6, 6, 5, 5, 5, 5, 5, 5, 4, 5, 6, 5, 5, 4, 4, 2, 5, 5, 5, 5, 6, 1, 4, 5, 5, 5, 4, 4, 5, 5, 5, 5, 5, 5, 6, 5, 4, 5, 6, 4, 5, 5, 5, 5, 5, 4, 4, 5, 2, 4, 5, 6, 4, 6, 5, 5, 5, 5, 4, 4, 4, 3, 5, 3, 1, 4, 4, 4, 4, 5, 5, 5, 4, 5, 4, 4, 5, 3, 5, 6, 5, 3, 4, 5, 5, 3, 3, 5, 5, 5, 5, 5, 5, 5, 5, 2, 5, 4, 4, 4, 5, 5, 4, 5, 5, 2, 5, 3, 2, 5, 4, 5, 4, 3, 4, 2, 4, 5, 5, 4, 4, 5, 5, 5, 4, 5, 5, 5, 4, 5, 5, 5, 5, 4, 4, 4, 6, 4, 5, 5, 2, 4, 3, 4, 5, 4, 5, 5, 5, 5, 4, 3, 5, 3, 5, 5, 5, 4, 5, 4, 4, 6, 4, 5, 3, 5, 5, 5, 5, 2, 5, 3, 4, 4, 3, 5, 6, 5, 4, 5, 4, 4, 4, 4, 5, 2, 3, 5, 5, 5, 5, 4, 3, 5, 5, 5, 6, 5, 4, 4, 5, 4, 5, 4, 5, 6, 5, 3, 3, 5, 5, 4, 4, 4, 3, 4, 6, 5, 5, 4, 5, 5, 5, 6, 5, 4, 6, 4, 4, 4, 4, 5, 5, 2, 4, 6, 6, 6, 4, 3, 3, 5, 4, 6, 3, 5, 5, 3, 3, 6, 5, 3, 4, 5, 4, 4, 4, 4, 4, 5, 5, 4, 5, 4, 6, 4, 4, 5, 5, 5, 4, 5, 2, 5, 2, 3, 2, 5, 5, 4, 5, 5, 5, 5, 5, 5, 5, 4, 5, 5, 6, 3, 4, 5, 4, 5, 4, 4, 3, 4, 5, 4, 4, 5, 3, 5, 6, 1, 5, 4, 5, 5, 3, 4, 4, 5, 2, 4,
Re: [R] Retrieve hypergeometric results in large scale
Hi Bert. This is not a homework. If I can do some basic programming in R like Perl, then I'll have a better chance to accomplish this task but the matrix concept is not quickly comprehensible... -- View this message in context: http://r.789695.n4.nabble.com/Retrieve-95-coverage-of-results-from-a-hypergeometric-distribution-tp4644683p4644703.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] Transform pairwise observations into a table
Hi, I have a table of pairs of individuals and a coefficient that belongs to the pair: ind1ind2coef 1 1 1 1 2 0.25 1 3 0.125 1 4 0.5 2 2 1 2 1 0.25 2 3 0.125 2 4 0.5 3 3 1 3 1 0.125 3 2 0.125 3 4 0.5 4 4 1 4 1 0.5 4 2 0.5 4 3 0.5 And I want to convert it to a matrix where each individual is both a row and a column and at the intersection of each pair is the coefficient that belongs to that pair: 1 2 3 4 1 1 0.250.1250.5 2 0.25 10.1250.5 3 0.1250.1251 0.5 4 0.50.5 0.5 1 If table() would allow me to specify something other than frequencies to fill the table with, it would be what I need. I tried a few different combinations of t() and unique() but none of it made enough sense to post as my starting code... I am just lost. Any help would be greatly appreciated. Thank you, AHJ -- View this message in context: http://r.789695.n4.nabble.com/Transform-pairwise-observations-into-a-table-tp4644706.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] Hmisc describe error
On Oct 1, 2012, at 9:33 AM, Bond, Stephen wrote: Describe fails for me with a message similar to what was an issue in 2008 and got fixed according to posts. R version 2.15.0 (2012-03-30) Copyright (C) 2012 The R Foundation for Statistical Computing ISBN 3-900051-07-0 Platform: i386-pc-mingw32/i386 (32-bit) # output truncated options(chmhelp = FALSE, help_type = text) .help.ESS - help options(STERM='iESS', editor='gnuclient.exe') load('prostate.sav') library(rms) Loading required package: Hmisc Loading required package: survival Loading required package: splines Hmisc library by Frank E Harrell Jr Type library(help='Hmisc'), ?Overview, or ?Hmisc.Overview') to see overall documentation. NOTE:Hmisc no longer redefines [.factor to drop unused levels when subsetting. To get the old behavior of Hmisc type dropUnusedLevels(). Attaching package: 'Hmisc' The following object(s) are masked from 'package:survival': untangle.specials The following object(s) are masked from 'package:base': format.pval, round.POSIXt, trunc.POSIXt, units Attaching package: 'rms' The following object(s) are masked from 'package:survival': Surv describe(prostate) Error in format(dates(x)) : could not find function dates Iget anerror when I use that code. Error in describe(prostate) : object 'prostate' not found And when I attempt loading it get: data(prostate) So you must have a different set of objects or packages loaded on your installation than I do. David Winsemius, MD Alameda, CA, USA __ 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] Transform pairwise observations into a table
Thank you. I had looked at xtabs but misunderstood the syntax. This is great. :) AHJ On Oct 1, 2012, at 12:53 PM, arun smartpink...@yahoo.com wrote: Hi, Try this: dat1-read.table(text= ind1 ind2 coef 1 1 1 1 2 0.25 1 3 0.125 1 4 0.5 2 2 1 2 1 0.25 2 3 0.125 2 4 0.5 3 3 1 3 1 0.125 3 2 0.125 3 4 0.5 4 4 1 4 1 0.5 4 2 0.5 4 3 0.5 ,sep=,header=TRUE) mat1-as.matrix(xtabs(coef~ind1+ind2,data=dat1)) #ind2 #ind1 1 2 3 4 # 1 1.000 0.250 0.125 0.500 #2 0.250 1.000 0.125 0.500 #3 0.125 0.125 1.000 0.500 #4 0.500 0.500 0.500 1.000 A.K. - Original Message - From: AHJ ahadjixenofon...@med.miami.edu To: r-help@r-project.org Cc: Sent: Monday, October 1, 2012 12:17 PM Subject: [R] Transform pairwise observations into a table Hi, I have a table of pairs of individuals and a coefficient that belongs to the pair: ind1ind2coef 111 120.25 130.125 140.5 221 210.25 230.125 240.5 331 310.125 320.125 340.5 441 410.5 420.5 430.5 And I want to convert it to a matrix where each individual is both a row and a column and at the intersection of each pair is the coefficient that belongs to that pair: 1 2 3 4 11 0.25 0.1250.5 20.25 1 0.1250.5 30.1250.12510.5 40.5 0.5 0.51 If table() would allow me to specify something other than frequencies to fill the table with, it would be what I need. I tried a few different combinations of t() and unique() but none of it made enough sense to post as my starting code... I am just lost. Any help would be greatly appreciated. Thank you, AHJ -- View this message in context: http://r.789695.n4.nabble.com/Transform-pairwise-observations-into-a-table-tp4644706.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] Retrieve hypergeometric results in large scale
If you have not already done so, stop what you are doing and work through the Introduction to R tutorial that ships with R (or other R tutorial on the web that you may prefer). The tutorials are written to help you climb the R learning curve much more efficiently than the fooling around that you appear to be doing now. -- Bert On Mon, Oct 1, 2012 at 8:31 AM, jas4710 wata...@post.com wrote: Hi Bert. This is not a homework. If I can do some basic programming in R like Perl, then I'll have a better chance to accomplish this task but the matrix concept is not quickly comprehensible... -- View this message in context: http://r.789695.n4.nabble.com/Retrieve-95-coverage-of-results-from-a-hypergeometric-distribution-tp4644683p4644703.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. -- 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] Hmisc describe error
The error is just a misleading error message. loading the data produced column sdate as str(prostate) $ sdate :Classes 'labelled', 'dates' atomic [1:502] 2778 2820 2933 2999 3002 ... .. ..- attr(*, format)= symbol ddmmmyy .. ..- attr(*, label)= chr Date on study prostate$sdate - as.Date(prostate$sdate,origin=1970-01-01) # fixes the problem and 'describe' runs just as a comment: maybe the R dataset should have sdate as a Date class. -- View this message in context: http://r.789695.n4.nabble.com/Hmisc-describe-error-tp4644708p4644711.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] Transform pairwise observations into a table
Hi, Try this: dat1-read.table(text= ind1 ind2 coef 1 1 1 1 2 0.25 1 3 0.125 1 4 0.5 2 2 1 2 1 0.25 2 3 0.125 2 4 0.5 3 3 1 3 1 0.125 3 2 0.125 3 4 0.5 4 4 1 4 1 0.5 4 2 0.5 4 3 0.5 ,sep=,header=TRUE) mat1-as.matrix(xtabs(coef~ind1+ind2,data=dat1)) # ind2 #ind1 1 2 3 4 # 1 1.000 0.250 0.125 0.500 #2 0.250 1.000 0.125 0.500 #3 0.125 0.125 1.000 0.500 #4 0.500 0.500 0.500 1.000 A.K. - Original Message - From: AHJ ahadjixenofon...@med.miami.edu To: r-help@r-project.org Cc: Sent: Monday, October 1, 2012 12:17 PM Subject: [R] Transform pairwise observations into a table Hi, I have a table of pairs of individuals and a coefficient that belongs to the pair: ind1 ind2 coef 1 1 1 1 2 0.25 1 3 0.125 1 4 0.5 2 2 1 2 1 0.25 2 3 0.125 2 4 0.5 3 3 1 3 1 0.125 3 2 0.125 3 4 0.5 4 4 1 4 1 0.5 4 2 0.5 4 3 0.5 And I want to convert it to a matrix where each individual is both a row and a column and at the intersection of each pair is the coefficient that belongs to that pair: 1 2 3 4 1 1 0.25 0.125 0.5 2 0.25 1 0.125 0.5 3 0.125 0.125 1 0.5 4 0.5 0.5 0.5 1 If table() would allow me to specify something other than frequencies to fill the table with, it would be what I need. I tried a few different combinations of t() and unique() but none of it made enough sense to post as my starting code... I am just lost. Any help would be greatly appreciated. Thank you, AHJ -- View this message in context: http://r.789695.n4.nabble.com/Transform-pairwise-observations-into-a-table-tp4644706.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] Retrieve hypergeometric results in large scale
Thanks Jeff~~~ In fact I do not know how to combine and extract vectors in R. ans-sort(dhyper(x, m, n, k),decreasing=TRUE) rbind(ans,cumsum(ans) will show the first point that exceeds 95% threshold. The problem is: *information is lost* I can no longer identify where are the first few elements from. e.g. for 10 numbers, maybe they are from 4,5,6,7 or for 100 numbers, from 45 to 68 So to append ID's to the data for later retrieval? rbind appears to do the job but not so exactly... -- View this message in context: http://r.789695.n4.nabble.com/Retrieve-95-coverage-of-results-from-a-hypergeometric-distribution-tp4644683p4644715.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] svyboxplot - library (survey)
try svyboxplot( dthage ~ factor( xspd2 ) , subset(nhis, mortstat==1) , col=gray80, varwidth=TRUE , ylab=Age at Death , xlab=SPD Status: 1-SPD, 2=No SPD ) On Mon, Oct 1, 2012 at 11:28 AM, Muhuri, Pradip (SAMHSA/CBHSQ) pradip.muh...@samhsa.hhs.gov wrote: Dear Anthony, ** ** Yes, I can follow the example code you have given. But, do you know from the code shown below (following Thomas Lumleys Complex Surveys) why I am getting the boxplot of dthage for just xspd=1, not xspd2=2? ** ** My intent is the make this code work so that I can generate similar plots on other continuous variable. ** ** Any help will be appreciated. ** ** Thanks, ** ** Pradip ** ** ** ** ** ** nhis - svydesign (id=~psu, strat=~stratum, weights=~wt8, data=tor, nest=TRUE) svyboxplot (dthage~xspd2, subset (nhis, mortstat==1), col=gray80, varwidth=TRUE, ylab=Age at Death, xlab=SPD Status: 1-SPD, 2=No SPD) ** ** Pradip K. Muhuri, PhD Statistician Substance Abuse Mental Health Services Administration The Center for Behavioral Health Statistics and Quality Division of Population Surveys 1 Choke Cherry Road, Room 2-1071 Rockville, MD 20857 Tel: 240-276-1070 Fax: 240-276-1260 e-mail: pradip.muh...@samhsa.hhs.gov The Center for Behavioral Health Statistics and Quality your feedback. Please click on the following link to complete a brief customer survey: http://cbhsqsurvey.samhsa.gov ** ** *From:* Anthony Damico [mailto:ajdam...@gmail.com] *Sent:* Monday, October 01, 2012 10:07 AM *To:* Muhuri, Pradip (SAMHSA/CBHSQ) *Cc:* R help *Subject:* Re: [R] svyboxplot - library (survey) ** ** using a slight modification of the example shown in ?svyboxplot # load survey library library(survey) # load example data data(api) # create an example svydesign dstrat - svydesign(id = ~1, strata = ~stype, weights = ~pw, data = apistrat, fpc = ~fpc) # set the plot window to display 1 plot x 2 plots par(mfrow=c(1,2)) # generate two example boxplots svyboxplot(enroll~stype,dstrat,all.outliers=TRUE) svyboxplot(enroll~1,dstrat) # done # alternative: not as nice # set the plot window to display 2 plots x 1 plot par(mfrow=c(2,1)) # generate two example boxplots svyboxplot(enroll~stype,dstrat,all.outliers=TRUE) svyboxplot(enroll~1,dstrat) # done On Mon, Oct 1, 2012 at 9:50 AM, Muhuri, Pradip (SAMHSA/CBHSQ) pradip.muh...@samhsa.hhs.gov wrote: Hello, I have used the library (survey) package for boxplots using the following code. Could anyone please tell me why I am getting only 1 boxplot instead of 2 boxplots (1-SPD, 2-No SPD). What changes in the following code would be required to get 2 boxplots in the same plot frame? Thanks, Pradip ### nhis - svydesign (id=~psu, strat=~stratum, weights=~wt8, data=tor, nest=TRUE) svyboxplot (dthage~xspd2, subset (nhis, mortstat==1), col=gray80, varwidth=TRUE, ylab=Age at Death, xlab=SPD Status: 1-SPD, 2=No SPD) Pradip K. Muhuri Statistician Substance Abuse Mental Health Services Administration The Center for Behavioral Health Statistics and Quality Division of Population Surveys 1 Choke Cherry Road, Room 2-1071 Rockville, MD 20857 Tel: 240-276-1070 Fax: 240-276-1260 e-mail: pradip.muh...@samhsa.hhs.gov The Center for Behavioral Health Statistics and Quality your feedback. Please click on the following link to complete a brief customer survey: http://cbhsqsurvey.samhsa.gov vide 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. ** ** [[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] Transform pairwise observations into a table
You might be amused by this alternative: with(dat1,matrix(coef[order(ind2,ind1)],nr=4,nc=4)) It works because matrices are vectors stored in column major order. Cheers, Bert On Mon, Oct 1, 2012 at 9:53 AM, arun smartpink...@yahoo.com wrote: Hi, Try this: dat1-read.table(text= ind1 ind2 coef 1 1 1 1 2 0.25 1 3 0.125 1 4 0.5 2 2 1 2 1 0.25 2 3 0.125 2 4 0.5 3 3 1 3 1 0.125 3 2 0.125 3 4 0.5 4 4 1 4 1 0.5 4 2 0.5 4 3 0.5 ,sep=,header=TRUE) mat1-as.matrix(xtabs(coef~ind1+ind2,data=dat1)) #ind2 #ind1 1 2 3 4 # 1 1.000 0.250 0.125 0.500 #2 0.250 1.000 0.125 0.500 #3 0.125 0.125 1.000 0.500 #4 0.500 0.500 0.500 1.000 A.K. - Original Message - From: AHJ ahadjixenofon...@med.miami.edu To: r-help@r-project.org Cc: Sent: Monday, October 1, 2012 12:17 PM Subject: [R] Transform pairwise observations into a table Hi, I have a table of pairs of individuals and a coefficient that belongs to the pair: ind1ind2coef 111 120.25 130.125 140.5 221 210.25 230.125 240.5 331 310.125 320.125 340.5 441 410.5 420.5 430.5 And I want to convert it to a matrix where each individual is both a row and a column and at the intersection of each pair is the coefficient that belongs to that pair: 1 2 3 4 11 0.25 0.1250.5 20.25 1 0.1250.5 30.1250.12510.5 40.5 0.5 0.51 If table() would allow me to specify something other than frequencies to fill the table with, it would be what I need. I tried a few different combinations of t() and unique() but none of it made enough sense to post as my starting code... I am just lost. Any help would be greatly appreciated. Thank you, AHJ -- View this message in context: http://r.789695.n4.nabble.com/Transform-pairwise-observations-into-a-table-tp4644706.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. -- 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] Hmisc describe error
On Oct 1, 2012, at 9:48 AM, stephenb wrote: The error is just a misleading error message. loading the data produced column sdate as str(prostate) $ sdate :Classes 'labelled', 'dates' atomic [1:502] 2778 2820 2933 2999 3002 ... .. ..- attr(*, format)= symbol ddmmmyy .. ..- attr(*, label)= chr Date on study prostate$sdate - as.Date(prostate$sdate,origin=1970-01-01) # fixes the problem and 'describe' runs just as a comment: maybe the R dataset should have sdate as a Date class. What R dataset? From what package? And it's already a dates-classed variable. -- David Winsemius, MD Alameda, CA, USA __ 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] Hmisc describe error
Just as a clarification: I downloaded 'prostate.sav' from F. Harrell's website. For some reason data(prostate) Warning message: In data(prostate) : data set 'prostate' not found I don't have the prostate data set as is. -Original Message- From: David Winsemius [mailto:dwinsem...@comcast.net] Sent: Monday, October 01, 2012 12:53 PM To: Bond, Stephen Cc: r-help@r-project.org Subject: Re: [R] Hmisc describe error On Oct 1, 2012, at 9:33 AM, Bond, Stephen wrote: Describe fails for me with a message similar to what was an issue in 2008 and got fixed according to posts. R version 2.15.0 (2012-03-30) Copyright (C) 2012 The R Foundation for Statistical Computing ISBN 3-900051-07-0 Platform: i386-pc-mingw32/i386 (32-bit) # output truncated options(chmhelp = FALSE, help_type = text) .help.ESS - help options(STERM='iESS', editor='gnuclient.exe') load('prostate.sav') library(rms) Loading required package: Hmisc Loading required package: survival Loading required package: splines Hmisc library by Frank E Harrell Jr Type library(help='Hmisc'), ?Overview, or ?Hmisc.Overview') to see overall documentation. NOTE:Hmisc no longer redefines [.factor to drop unused levels when subsetting. To get the old behavior of Hmisc type dropUnusedLevels(). Attaching package: 'Hmisc' The following object(s) are masked from 'package:survival': untangle.specials The following object(s) are masked from 'package:base': format.pval, round.POSIXt, trunc.POSIXt, units Attaching package: 'rms' The following object(s) are masked from 'package:survival': Surv describe(prostate) Error in format(dates(x)) : could not find function dates Iget anerror when I use that code. Error in describe(prostate) : object 'prostate' not found And when I attempt loading it get: data(prostate) So you must have a different set of objects or packages loaded on your installation than I do. David Winsemius, MD Alameda, CA, USA __ 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] Hmisc describe error
On Oct 1, 2012, at 10:27 AM, David Winsemius wrote: On Oct 1, 2012, at 9:48 AM, stephenb wrote: The error is just a misleading error message. loading the data produced column sdate as str(prostate) $ sdate :Classes 'labelled', 'dates' atomic [1:502] 2778 2820 2933 2999 3002 ... .. ..- attr(*, format)= symbol ddmmmyy .. ..- attr(*, label)= chr Date on study prostate$sdate - as.Date(prostate$sdate,origin=1970-01-01) # fixes the problem and 'describe' runs just as a comment: maybe the R dataset should have sdate as a Date class. What R dataset? From what package? And it's already a dates-classed variable. Not only that, but there is a message at the top of Harrell's Datasets page that says: For R users of the prostate dataset, put library(chron) into effect to handle date variables. A simpler approach is to just convert the one date variable to the built-in R format by running the command prostate$sdate - as.Date(prostate$sdate). -- David Winsemius, MD Alameda, CA, USA __ 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] Transform pairwise observations into a table
Hello, Try the following. dat - read.table(text= ind1 ind2 coef 1 1 1 1 2 0.25 1 3 0.125 1 4 0.5 2 2 1 2 1 0.25 2 3 0.125 2 4 0.5 3 3 1 3 1 0.125 3 2 0.125 3 4 0.5 4 4 1 4 1 0.5 4 2 0.5 4 3 0.5 , header=TRUE) dat reshape(dat, v.names = coef, idvar = ind1, timevar = ind2, direction = wide) Hope this helps, Rui Barradas Em 01-10-2012 17:17, AHJ escreveu: Hi, I have a table of pairs of individuals and a coefficient that belongs to the pair: ind1ind2coef 1 1 1 1 2 0.25 1 3 0.125 1 4 0.5 2 2 1 2 1 0.25 2 3 0.125 2 4 0.5 3 3 1 3 1 0.125 3 2 0.125 3 4 0.5 4 4 1 4 1 0.5 4 2 0.5 4 3 0.5 And I want to convert it to a matrix where each individual is both a row and a column and at the intersection of each pair is the coefficient that belongs to that pair: 1 2 3 4 1 1 0.250.1250.5 2 0.25 10.1250.5 3 0.1250.1251 0.5 4 0.50.5 0.5 1 If table() would allow me to specify something other than frequencies to fill the table with, it would be what I need. I tried a few different combinations of t() and unique() but none of it made enough sense to post as my starting code... I am just lost. Any help would be greatly appreciated. Thank you, AHJ -- View this message in context: http://r.789695.n4.nabble.com/Transform-pairwise-observations-into-a-table-tp4644706.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] Retrieve hypergeometric results in large scale
Hello, See the differences. k - 3 p - 0.95 m - 90; n - 10 dhyper(0:k, m, n, k) # Prob(X = x), with x = 0:k phyper(0:k, m, n, k) # Prob(X = x) # quantiles, what you want qhyper(p, m, n, k) # inverse of phyper m - 50; n - 50 dhyper(0:k, m, n, k) phyper(0:k, m, n, k) qhyper(p, m, n, k) In your original post (op) you gave the output of dhyper. I would find phyper much better, it shows that 0.95 is between the last two values of phyper. Apparently you want its inverse, qhyper, but are assuming there are such things as decimals. The hypergeometric is a discrete distribution, so, from the help page, The quantile is defined as the smallest value /x/ such that /F(x) ? p/, where /F/ is the distribution function. And, R is surely more competent at distribution functions than Perl. Stick to it and read An Introduction to R, file R-intro.pdf in the doc directory of your R installation, Chapter 8. All computer languages have some learning time and as far as statistics is concerned, learning R pays. Hope this helps, Rui Barradas Em 01-10-2012 16:28, jas4710 escreveu: Thanks Jeff The documentation pages, if I haven't missed any crucial points, illustrate how to get probability and cumulative probability values. I can first retrieve the data structures and use Perl (I don't know how to use R...) to sort the derived ratios and sum the probability values until the cumulative probability exceeds 95%. Well, I just don't know whether such seemingly routine procedures have already been implemented... Thanks again! -- View this message in context: http://r.789695.n4.nabble.com/Retrieve-95-coverage-of-results-from-a-hypergeometric-distribution-tp4644683p4644701.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.
[R] ffbase, help with %in%
Hello to everyone. I'm trying to use the %in% to match to vectors in ff format. a-as.ff(data[,1]) %in% fire$fecha aff (open) logical length=3653 (3653) [1][2][3][4][5][6][7][8][3646] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE : FALSE [3647] [3648] [3649] [3650] [3651] [3652] [3653] FALSE FALSE FALSE FALSE FALSE FALSE FALSE Here you see a part of the data: data[1:20,] (just a sample, data has 3653 obs) fecha juliano altura UTM.E UTM.N 1 1990-07-01 182 15 248500 6239500 2 1990-07-02 183 15 248500 6239500 3 1990-07-03 184 15 248500 6239500 4 1990-07-04 185 15 248500 6239500 5 1990-07-05 186 15 248500 6239500 6 1990-07-06 187 15 248500 6239500 7 1990-07-07 188 15 248500 6239500 8 1990-07-08 189 15 248500 6239500 9 1990-07-09 190 15 248500 6239500 10 1990-07-10 191 15 248500 6239500 11 1990-07-11 192 15 248500 6239500 12 1990-07-12 193 15 248500 6239500 13 1990-07-13 194 15 248500 6239500 14 1990-07-14 195 15 248500 6239500 15 1990-07-15 196 15 248500 6239500 16 1990-07-16 197 15 248500 6239500 17 1990-07-17 198 15 248500 6239500 18 1990-07-18 199 15 248500 6239500 19 1990-07-19 200 15 248500 6239500 20 1990-07-20 201 15 248500 6239500 fire$fecha[1:20,] [1] 1984-11-08 1984-11-08 1984-11-09 1984-11-09 1984-11-09 [6] 1984-11-10 1984-11-10 1984-11-11 1984-11-11 1984-11-11 [11] 1984-11-11 1984-11-11 1984-11-11 1984-11-12 1984-11-12 [16] 1984-11-13 1984-11-13 1984-11-13 1984-11-14 1984-11-14 to see if a got any match: table.ff(a) FALSE TRUE 1687 1966 Mensajes de aviso perdidosIn if (useNA == no) c(NA, NaN) : la condición tiene longitud 1 y sólo el primer elemento será usado in a regular data.frame I use data[a,] to extract the rows that a == TRUE, but when i do this in a ffdf i get this error: data[a,]Error: vmode(index) == integer is not TRUE I'm just learning how to use the ff package so, obviously I'm missing something If any of you knows how to solve this, please teach me. Thank you so much. Lucas. [[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] enquiry
Hello, Or maybe %in%. x - Sys.Date() + 1:10 y - Sys.Date() + 7:15 x %in% y # logical index into 'x' x[x %in% y] # common x[!x %in% y] # not common setdiff(x, y) # doesn't keep class Date Hope this helps, Rui Barradas Em 01-10-2012 15:57, R. Michael Weylandt escreveu: On Monday, October 1, 2012, Karan Anand wrote: hi, I am new to using R ,I have 2 datasets with dates common in them ,how can i take out the common dates within them. Depending on what you mean by 'out' either merge() or setdiff() RMW karan [[alternative HTML version deleted]] __ R-help@r-project.org javascript:; 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-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] modified Ward method
You asked this question before. I think you need to tell us what modified ward method you are talking about and what you are trying to accomplish. I did a quick Google Scholar search and turned up only a few papers using that phrase, but they were not defining it in the same way. -- David L Carlson Associate Professor of Anthropology Texas AM University College Station, TX 77843-4352 -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-bounces@r- project.org] On Behalf Of eliza botto Sent: Saturday, September 29, 2012 9:29 PM To: r-help@r-project.org Subject: [R] modified Ward method Dear useRs,Is there a command or package in R to do heirarchical clustering by modified ward method??regardseliza [[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] calculating probability from the density function
You may find it easier to use the logspline density fits (logspline package) rather than the kernel density estimators for this. On Mon, Oct 1, 2012 at 7:46 AM, Eugene Kanshin kanshin...@gmail.com wrote: Hello, I have a data x with normal (or very close to normal) distribution, I can plot a density distribution with density(x,...). My question is is there any way to calculate an area under this distribution (=probability) for particular range of x values, let's say for x from 0 to 2? I was not able to find any kind of simple procedure to do this. Thanks in advance for your help, Evgeny. [[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.
[R] (senza oggetto)
hi, can I help me to cancel my name in this list? 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] (senza oggetto)
On Oct 1, 2012, at 20:28 , MYRIAM TABASSO wrote: hi, can I help me to cancel my name in this list? Sure, just go to the mailman site listed in the footer (look near bottom of page). 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. -- 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] ffbase, help with %in%
Hi Lucas I don't know the ff package very well but here is what I found. Maybe there is a clue in here z - as.Date(1970-01-01)+1:10 zff - as.ff(z) z %in% zff [1] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE #but z %in% zff[1:length(zff),] [1] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE # and also z %in% zff[1:10] [1] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE That seems to be because an ff object is a list str(zff) list() - attr(*, physical)=Class 'ff_pointer' externalptr ..- attr(*, vmode)= chr double ..- attr(*, maxlength)= int 10 ..- attr(*, pattern)= chr clone ..- attr(*, filename)= chr /private/var/folders/ub/ubvWLUkKHf8WAywv5rmtcE+++TI/-Tmp-/RtmpmMKHRx/clone406828b135cc.ff ..- attr(*, pagesize)= int 65536 ..- attr(*, finalizer)= chr close ..- attr(*, finonexit)= logi TRUE ..- attr(*, readonly)= logi FALSE ..- attr(*, caching)= chr mmnoflush - attr(*, virtual)= list() ..- attr(*, Length)= int 10 ..- attr(*, Symmetric)= logi FALSE ..- attr(*, ramclass)= chr Date - attr(*, class) = chr [1:2] ff_vector ff If you make a data frame e.g df - data.frame(date=z,letter=letters[1:10]) # this doesn't work as.ff(df[1:10,1])%in% z logical(0) # but this does as.ff(df[1:10,1]%in% z) ff (open) logical length=10 (10) [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE The question is: In your example, do you want a to be a ff object or the first column of data to be an ff object (I could get this to work - propably for the reasons above : df - data.frame(date=as.ff(z),letter=letters[1:10]) Error in as.data.frame.default(x[[i]], optional = TRUE, stringsAsFactors = stringsAsFactors) : cannot coerce class 'c(ff_vector, ff)' into a data.frame On 1 October 2012 20:07, Lucas Chaparro lpchaparro...@gmail.com wrote: Hello to everyone. I'm trying to use the %in% to match to vectors in ff format. a-as.ff(data[,1]) %in% fire$fecha aff (open) logical length=3653 (3653) [1][2][3][4][5][6][7][8][3646] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE : FALSE [3647] [3648] [3649] [3650] [3651] [3652] [3653] FALSE FALSE FALSE FALSE FALSE FALSE FALSE Here you see a part of the data: data[1:20,] (just a sample, data has 3653 obs) fecha juliano altura UTM.E UTM.N 1 1990-07-01 182 15 248500 6239500 2 1990-07-02 183 15 248500 6239500 3 1990-07-03 184 15 248500 6239500 4 1990-07-04 185 15 248500 6239500 5 1990-07-05 186 15 248500 6239500 6 1990-07-06 187 15 248500 6239500 7 1990-07-07 188 15 248500 6239500 8 1990-07-08 189 15 248500 6239500 9 1990-07-09 190 15 248500 6239500 10 1990-07-10 191 15 248500 6239500 11 1990-07-11 192 15 248500 6239500 12 1990-07-12 193 15 248500 6239500 13 1990-07-13 194 15 248500 6239500 14 1990-07-14 195 15 248500 6239500 15 1990-07-15 196 15 248500 6239500 16 1990-07-16 197 15 248500 6239500 17 1990-07-17 198 15 248500 6239500 18 1990-07-18 199 15 248500 6239500 19 1990-07-19 200 15 248500 6239500 20 1990-07-20 201 15 248500 6239500 fire$fecha[1:20,] [1] 1984-11-08 1984-11-08 1984-11-09 1984-11-09 1984-11-09 [6] 1984-11-10 1984-11-10 1984-11-11 1984-11-11 1984-11-11 [11] 1984-11-11 1984-11-11 1984-11-11 1984-11-12 1984-11-12 [16] 1984-11-13 1984-11-13 1984-11-13 1984-11-14 1984-11-14 to see if a got any match: table.ff(a) FALSE TRUE 1687 1966 Mensajes de aviso perdidosIn if (useNA == no) c(NA, NaN) : la condición tiene longitud 1 y sólo el primer elemento será usado in a regular data.frame I use data[a,] to extract the rows that a == TRUE, but when i do this in a ffdf i get this error: data[a,]Error: vmode(index) == integer is not TRUE I'm just learning how to use the ff package so, obviously I'm missing something If any of you knows how to solve this, please teach me. Thank you so much. Lucas. [[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. -- Christiaan Pauw Nova Institute www.nova.org.za __ 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 nls regression fit
You have not specified a nonlinear formula. There are no parameters to estimate in the formula you provide, y1~dist. What is the nonlinear relation you are trying to fit? Look at the help file for nls to see some examples worked. ?nls Jean Gyanendra Pokharel gyanendra.pokha...@gmail.com wrote on 10/01/2012 10:27:23 AM: Hi all, I got following problem in fitting the data. Any kind of suggestions are welcome beta - 3.5 d - seq(0.1,62.5,0.1) y - exp(-beta*d) y1 - y x - read.table(epidist.txt, header = TRUE) data.nls - as.data.frame(cbind(y1,x)) #attach(data.nls) nls.fit - nls(y1~dist,data.nls) Error in cll[[1L]] : object of type 'symbol' is not subsettable Best [[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] clustering spline-based models
Hello playeRs! I'm working on a project for a client. She's modeling hormone levels periodically, and trying to develop a model and fit her data to that model, and subsequently she's trying to cluster individuals based on how well each fits the model. I've been looking at grofit for this, but it doesn't appear that it can do the sort of post-hoc analysis I'm looking for, although it can use a model-free spline fit nicely (reportedly). Does anyone have any packages that can accomplish what I'm looking to accomplish? Thanks in advance for any help you can offer, Wyatt [[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] Retrieve hypergeometric results in large scale
order() is usually a lot more useful than sort(), since, as you noticed, sort() drops information about where each element in its output came from. Your example was incomplete so I made up one which I think is similar. n - 10 ; p - 0.7 ; k - 0:n ; d - dbinom(k, n, p) plot(k, d) # density of binomial over its domain If you want the indices of the largest density values whose cumulative sum is less than 0.95 you ord - order(d, decreasing=TRUE) # indices such that d[ord] is in decreasing order big - ord[cumsum(d[ord]) 0.95] data.frame(big, d=d[big], cumsum=cumsum(d[big])) big dcumsum 1 8 0.2668279 0.2668279 2 9 0.2334744 0.5003024 3 7 0.2001209 0.7004233 4 10 0.1210608 0.8214841 5 6 0.1029193 0.9244035 points(cex=2, k[big], d[big]) If you want to include the index of the density value that puts you just over 0.95 first find the complement of the desired indices and use setdiff to compute its complement. E.g., ord - order(d) little - ord[cumsum(d[ord]) 0.05] big - setdiff(seq_along(d), little) # difference of two sets of numbers big [1] 5 6 7 8 9 10 Bill Dunlap Spotfire, TIBCO Software wdunlap tibco.com -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of jas4710 Sent: Monday, October 01, 2012 9:59 AM To: r-help@r-project.org Subject: Re: [R] Retrieve hypergeometric results in large scale Thanks Jeff~~~ In fact I do not know how to combine and extract vectors in R. ans-sort(dhyper(x, m, n, k),decreasing=TRUE) rbind(ans,cumsum(ans) will show the first point that exceeds 95% threshold. The problem is: *information is lost* I can no longer identify where are the first few elements from. e.g. for 10 numbers, maybe they are from 4,5,6,7 or for 100 numbers, from 45 to 68 So to append ID's to the data for later retrieval? rbind appears to do the job but not so exactly... -- View this message in context: http://r.789695.n4.nabble.com/Retrieve-95-coverage-of- results-from-a-hypergeometric-distribution-tp4644683p4644715.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] svyboxplot - library (survey)
The documentation says The grouping variable in svyboxplot, if present, must be a factor -thomas On Tue, Oct 2, 2012 at 4:28 AM, Muhuri, Pradip (SAMHSA/CBHSQ) pradip.muh...@samhsa.hhs.gov wrote: Dear Anthony, Yes, I can follow the example code you have given. But, do you know from the code shown below (following Thomas Lumley's Complex Surveys) why I am getting the boxplot of dthage for just xspd=1, not xspd2=2? My intent is the make this code work so that I can generate similar plots on other continuous variable. Any help will be appreciated. Thanks, Pradip nhis - svydesign (id=~psu, strat=~stratum, weights=~wt8, data=tor, nest=TRUE) svyboxplot (dthage~xspd2, subset (nhis, mortstat==1), col=gray80, varwidth=TRUE, ylab=Age at Death, xlab=SPD Status: 1-SPD, 2=No SPD) Pradip K. Muhuri, PhD Statistician Substance Abuse Mental Health Services Administration The Center for Behavioral Health Statistics and Quality Division of Population Surveys 1 Choke Cherry Road, Room 2-1071 Rockville, MD 20857 Tel: 240-276-1070 Fax: 240-276-1260 e-mail: pradip.muh...@samhsa.hhs.govmailto:pradip.muh...@samhsa.hhs.gov The Center for Behavioral Health Statistics and Quality your feedback. Please click on the following link to complete a brief customer survey: http://cbhsqsurvey.samhsa.govhttp://cbhsqsurvey.samhsa.gov/ From: Anthony Damico [mailto:ajdam...@gmail.com] Sent: Monday, October 01, 2012 10:07 AM To: Muhuri, Pradip (SAMHSA/CBHSQ) Cc: R help Subject: Re: [R] svyboxplot - library (survey) using a slight modification of the example shown in ?svyboxplot # load survey library library(survey) # load example data data(api) # create an example svydesign dstrat - svydesign(id = ~1, strata = ~stype, weights = ~pw, data = apistrat, fpc = ~fpc) # set the plot window to display 1 plot x 2 plots par(mfrow=c(1,2)) # generate two example boxplots svyboxplot(enroll~stype,dstrat,all.outliers=TRUE) svyboxplot(enroll~1,dstrat) # done # alternative: not as nice # set the plot window to display 2 plots x 1 plot par(mfrow=c(2,1)) # generate two example boxplots svyboxplot(enroll~stype,dstrat,all.outliers=TRUE) svyboxplot(enroll~1,dstrat) # done On Mon, Oct 1, 2012 at 9:50 AM, Muhuri, Pradip (SAMHSA/CBHSQ) pradip.muh...@samhsa.hhs.govmailto:pradip.muh...@samhsa.hhs.gov wrote: Hello, I have used the library (survey) package for boxplots using the following code. Could anyone please tell me why I am getting only 1 boxplot instead of 2 boxplots (1-SPD, 2-No SPD). What changes in the following code would be required to get 2 boxplots in the same plot frame? Thanks, Pradip ### nhis - svydesign (id=~psu, strat=~stratum, weights=~wt8, data=tor, nest=TRUE) svyboxplot (dthage~xspd2, subset (nhis, mortstat==1), col=gray80, varwidth=TRUE, ylab=Age at Death, xlab=SPD Status: 1-SPD, 2=No SPD) Pradip K. Muhuri Statistician Substance Abuse Mental Health Services Administration The Center for Behavioral Health Statistics and Quality Division of Population Surveys 1 Choke Cherry Road, Room 2-1071 Rockville, MD 20857 Tel: 240-276-1070 Fax: 240-276-1260 e-mail: pradip.muh...@samhsa.hhs.govmailto:pradip.muh...@samhsa.hhs.gov The Center for Behavioral Health Statistics and Quality your feedback. Please click on the following link to complete a brief customer survey: http://cbhsqsurvey.samhsa.gov vide commented, minimal, self-contained, reproducible code. __ R-help@r-project.orgmailto: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. -- 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.
Re: [R] Error messages when attempting to calculate polychoric correlation matrices
Dear Kevin, From ?polychor (in the polycor package): polychor(x, y, ML = FALSE, control = list(), std.err = FALSE, maxcor=.) Arguments x a contingency table of counts or an ordered categorical variable; the latter can be numeric, logical, a factor, or an ordered factor, but if a factor, its levels should be in proper order. y if x is a variable, a second ordered categorical variable. So polychor doesn't take a data frame as its argument, and computes only one polychoric correlation at a time. (It generally helps in R, and even more generally, to read the documentation.) Thus, e.g., with(data.matrix18, polychor(SABAS02, SABAS06)) [1] 0.1811532 with(data.matrix18, polychor(SABAS02, SABAS06, ML=TRUE)) [1] 0.1817071 If you want to compute a matrix of polychoric correlations using the polycor package, then you can use the hetcor function -- see ?hetcor -- but you'll then have to change your variables from numeric to factors or ordered factors, or hetcor will compute Pearson pairwise correlations among them. You could use, for example for (j in 1:ncol(data.matrix18)) data.matrix18[, j] - as.factor(data.matrix18[, j]) I hope this helps, John --- John Fox Senator McMaster Professor of Social Statistics Department of Sociology McMaster University Hamilton, Ontario, Canada -Original Message- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of Kevin Cheung Sent: Monday, October 01, 2012 12:19 PM To: 'r-help@r-project.org' Subject: [R] Error messages when attempting to calculate polychoric correlation matrices Dear R users, I am a psychology postgraduate student who is relatively new to using R. I am currently developing a psychometric scale and have run into a few problems when using R to calculate a polychoric correlation matrix for my dataset. I am trying to produce a polychoric correlation matrix for calculating ordinal reliability estimates (eg. Alpha, omega).The set consists of 439 observations of 18 variables. Each of the 18 variables was measured using a 6 point Likert scale but some variables only feature responses across 5 ordinal steps because there were no observed responses for one extreme of the item. The dataset has no missing data following multiple imputation and preliminary analyses were conducted using SPSS. I have included a command that will reproduce the dataset below, apologies for the size of the command: data.matrix18 - structure(list(SABAS02 = structure(c(6, 2, 5, 6, 4, 6, 5, 3, 5, 1, 4, 4, 6, 5, 5, 5, 4, 4, 4, 3, 4, 6, 5, 2, 5, 4, 5, 4, 4, 5, 5, 5, 2, 5, 4, 5, 5, 5, 6, 3, 6, 4, 5, 5, 5, 5, 4, 4, 5, 5, 5, 4, 5, 5, 4, 6, 6, 5, 6, 4, 4, 5, 4, 5, 5, 5, 6, 4, 4, 5, 4, 6, 5, 6, 5, 6, 6, 5, 4, 6, 6, 5, 3, 4, 2, 4, 5, 5, 6, 2, 6, 5, 4, 6, 2, 5, 3, 5, 4, 4, 5, 4, 4, 3, 4, 4, 6, 5, 6, 6, 4, 4, 5, 5, 5, 6, 5, 5, 5, 5, 4, 5, 5, 5, 4, 4, 5, 6, 4, 5, 5, 4, 5, 5, 4, 5, 5, 4, 5, 5, 5, 4, 2, 5, 4, 5, 5, 5, 5, 5, 5, 4, 5, 4, 4, 5, 4, 6, 2, 4, 6, 4, 4, 5, 4, 4, 6, 6, 5, 3, 3, 2, 4, 5, 5, 5, 5, 6, 2, 5, 5, 4, 6, 4, 5, 6, 5, 4, 5, 5, 4, 4, 5, 4, 5, 5, 3, 5, 4, 3, 5, 4, 5, 5, 5, 5, 5, 5, 5, 5, 4, 4, 5, 4, 4, 5, 5, 6, 4, 3, 5, 6, 5, 4, 5, 5, 2, 5, 5, 4, 2, 4, 5, 5, 5, 5, 5, 3, 4, 5, 2, 4, 3, 4, 5, 5, 5, 5, 5, 5, 5, 6, 6, 5, 4, 5, 6, 4, 5, 3, 5, 6, 5, 4, 5, 5, 5, 5, 5, 4, 4, 5, 6, 5, 4, 5, 5, 6, 6, 6, 4, 3, 4, 4, 5, 6, 4, 5, 5, 4, 6, 5, 5, 4, 5, 4, 4, 4, 4, 5, 3, 5, 5, 4, 3, 4, 5, 4, 4, 5, 4, 5, 5, 3, 5, 6, 5, 5, 5, 5, 5, 4, 3, 5, 4, 4, 5, 5, 6, 4, 5, 4, 2, 6, 5, 5, 5, 5, 5, 5, 5, 6, 6, 5, 5, 5, 4, 5, 5, 6, 5, 4, 4, 5, 4, 4, 5, 5, 4, 6, 4, 4, 5, 6, 5, 5, 4, 5, 4, 5, 4, 6, 5, 5, 5, 4, 5, 6, 5, 5, 5, 5, 5, 4, 3, 4, 5, 4, 4, 6, 4, 4, 4, 2, 5, 4, 4, 6, 2, 6, 4, 5, 5, 6, 5, 6, 4, 4, 3, 4, 5, 5, 6, 5, 5, 5, 5, 4, 5, 3, 5, 5, 5, 5, 6, 3, 5, 6, 5, 6, 4, 5, 4, 4, 5, 6, 2, 5, 6), value.labels = structure(c(6, 5, 4, 3, 2, 1), .Names = c(Strongly Agree, Agree, Slightly Agree, Slightly Disagree, Disagree, Strongly Disagree))), SABAS06 = structure(c(6, 6, 6, 6, 6, 6, 5, 5, 6, 3, 6, 6, 6, 5, 5, 6, 5, 6, 6, 5, 5, 6, 6, 6, 6, 5, 6, 3, 6, 5, 5, 5, 5, 6, 6, 6, 3, 6, 6, 5, 5, 4, 5, 6, 6, 5, 6, 5, 6, 6, 6, 6, 6, 6, 5, 5, 6, 6, 6, 6, 6, 6, 5, 6, 6, 6, 6, 6, 4, 6, 5, 5, 6, 6, 6, 6, 6, 6, 3, 5, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 5, 5, 6, 6, 6, 5, 5, 5, 6, 5, 6, 5, 6, 6, 6, 6, 5, 5, 6, 5, 5, 6, 2, 5, 5, 6, 5, 5, 5, 6, 5, 6, 5, 5, 6, 6, 6, 5, 5, 5, 6, 6, 5, 6, 6, 6, 4, 6, 6, 5, 4, 6, 5, 5, 6, 5, 5, 5, 6, 6, 6, 6, 5, 5, 6, 6, 5, 6, 4, 6, 6, 6, 6, 5, 5, 6, 4, 5, 6, 5, 6, 5, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 5, 5, 6, 6, 5, 5, 6, 6, 5, 6, 5, 6, 5, 4, 6, 5, 6, 6, 5, 6, 5, 6, 5, 5, 5, 5, 6, 6, 5, 6, 5, 4, 6, 5, 6, 6, 5, 5, 6, 5, 5, 6, 6, 6, 6, 5, 5, 6, 6, 5, 5, 6, 6, 4, 5, 5, 6, 6, 6, 6, 4, 6, 6, 6, 6, 6, 5, 6, 6, 6, 6, 5, 6, 6, 6, 5, 5, 6, 5, 6, 6, 6, 5, 2, 6, 6, 5, 5, 6, 5, 6, 6, 5, 6, 5, 6, 5, 6, 6, 6, 5, 5, 4, 5, 6, 6, 5, 5, 4, 4, 6, 6,
Re: [R] Transform pairwise observations into a table
HI AHJ, No problem. One more way in addition to reshape() (Rui's suggestion) to get the same result. library(reshape) as.matrix(cast(melt(dat1,id=c(ind1,ind2)),ind1~ind2,value=value)) # 1 2 3 4 #1 1.000 0.250 0.125 0.5 #2 0.250 1.000 0.125 0.5 #3 0.125 0.125 1.000 0.5 #4 0.500 0.500 0.500 1.0 A.K. - Original Message - From: Hadjixenofontos, Athena ahadjixenofon...@med.miami.edu To: arun smartpink...@yahoo.com Cc: R help r-help@r-project.org Sent: Monday, October 1, 2012 12:59 PM Subject: Re: [R] Transform pairwise observations into a table Thank you. I had looked at xtabs but misunderstood the syntax. This is great. :) AHJ On Oct 1, 2012, at 12:53 PM, arun smartpink...@yahoo.com wrote: Hi, Try this: dat1-read.table(text= ind1 ind2 coef 1 1 1 1 2 0.25 1 3 0.125 1 4 0.5 2 2 1 2 1 0.25 2 3 0.125 2 4 0.5 3 3 1 3 1 0.125 3 2 0.125 3 4 0.5 4 4 1 4 1 0.5 4 2 0.5 4 3 0.5 ,sep=,header=TRUE) mat1-as.matrix(xtabs(coef~ind1+ind2,data=dat1)) # ind2 #ind1 1 2 3 4 # 1 1.000 0.250 0.125 0.500 #2 0.250 1.000 0.125 0.500 #3 0.125 0.125 1.000 0.500 #4 0.500 0.500 0.500 1.000 A.K. - Original Message - From: AHJ ahadjixenofon...@med.miami.edu To: r-help@r-project.org Cc: Sent: Monday, October 1, 2012 12:17 PM Subject: [R] Transform pairwise observations into a table Hi, I have a table of pairs of individuals and a coefficient that belongs to the pair: ind1 ind2 coef 1 1 1 1 2 0.25 1 3 0.125 1 4 0.5 2 2 1 2 1 0.25 2 3 0.125 2 4 0.5 3 3 1 3 1 0.125 3 2 0.125 3 4 0.5 4 4 1 4 1 0.5 4 2 0.5 4 3 0.5 And I want to convert it to a matrix where each individual is both a row and a column and at the intersection of each pair is the coefficient that belongs to that pair: 1 2 3 4 1 1 0.25 0.125 0.5 2 0.25 1 0.125 0.5 3 0.125 0.125 1 0.5 4 0.5 0.5 0.5 1 If table() would allow me to specify something other than frequencies to fill the table with, it would be what I need. I tried a few different combinations of t() and unique() but none of it made enough sense to post as my starting code... I am just lost. Any help would be greatly appreciated. Thank you, AHJ -- View this message in context: http://r.789695.n4.nabble.com/Transform-pairwise-observations-into-a-table-tp4644706.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] Reading labels for very large heatmaps
Hello, I have a large (919 X 919), hierarchically clustered heatmap that I would like to read the labels off of. I have tried saving the figure in pdf format and enlarging it until I can see the labels, but if I make the labels small enough to read (so that they don't overlap) using cexRow and cexCol, they do not appear in the pdf. The limit seems to be anything below cexRow or Col = 0.06. Is there a way to have smaller labels visible in the pdf? Is this an issue with resolution? I could get by with just a quarter of the heatmap (i.e. half of a row and half of a column) so that the labels don't have to be so small, but the heatmap must be clustered before this is done. I tried to cut the dendrogram at just the halfway point of the heatmap, but was not successful. Any suggestions? Thanks, Jacob -- View this message in context: http://r.789695.n4.nabble.com/Reading-labels-for-very-large-heatmaps-tp4644739.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] (no subject)
Hello, I am a new R -user and request your help for the following problem. I need to merge two dataset of longitudinal study which has two column (id and respose) common. when I used merge option to join the datas side be side, because of the repeated subject id, I got larger data set which is not accurate. I would like to connect twi data sets by id and response in such a way that data are connected by same id and response type and if the same subject has less data point in one data set, then produce NA. Sample data sets is attached. Bibek __ 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] mlogit and model-based recursive partitioning
Hello: Has anyone tried to model-based recursive partition (using mob from package party; thanks Achim and colleagues) a data set based on a multinomial logit model (using mlogit from package mlogit; thanks Yves)? I attempted to do so, but there are at least two reasons why I could not. First, in mob I am not quite sure that a model of class StatModel exists for mlogit models. Second, as mlogit uses the pipe character | to specify the model, I wonder how this would interact with mob which uses pipe to differentiate between explanatory and segmentation variables. An example (not working) of what I would like to accomplish follows below. Thanks a lot. Tudor library(party) library(mlogit) data(Fishing, package = mlogit) Fish - mlogit.data(Fishing, varying = c(2:9), shape = wide, choice = mode) # FIT AN mlogit MODEL m1 - mlogit(mode ~ price + catch, data=Fish) # THE DESIRED END RESULT: SEGMENT m1 BASED ON INCOME AND/OR OTHER POSSIBLE COVARIATES -- View this message in context: http://r.789695.n4.nabble.com/mlogit-and-model-based-recursive-partitioning-tp4644743.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] glmmPQL and spatial autocorrelation
Hi all, I am analyzing data on habitat utilization of seals in the Southern Ocean. My data show spatial autocorrelation, which I'm interested in incorporating into my model. I am trying to model the presence of dives (versus simulated pseudo-absences) using a binomial generalized binomial model (glmmPQL), since I can incorporate the autocorrelation structure to the model using that package. However, I'm trying to come up with ways to evaluate which semivariogram is the most adequate for my model, but not having AICs complicated things a bit for me. Any suggestions on how to best achieve this? Thanks in advance Luis A. Huckstadt, Ph.D. Department of Ecology and Evolutionary Biology University of California Santa Cruz Long Marine Lab 100 Shaffer Road Santa Cruz, CA 95060 [[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)
It will be great if you can reproduce some code and data you are working on. Bhupendrasinh Thakre On Oct 1, 2012, at 4:25 PM, bibek sharma mbhpat...@gmail.com wrote: Hello, I am a new R -user and request your help for the following problem. I need to merge two dataset of longitudinal study which has two column (id and respose) common. when I used merge option to join the datas side be side, because of the repeated subject id, I got larger data set which is not accurate. I would like to connect twi data sets by id and response in such a way that data are connected by same id and response type and if the same subject has less data point in one data set, then produce NA. Sample data sets is attached. Bibek __ 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] Ifelse Execution
Hi Everyone, I am trying to run a time based query and need some of your help. Not much of data or packages. Just a simple one. Query I am trying to execute. ifelse ((as.numeric(as.POSIXct(2012-10-01 20:38:00))), (rnorm(1,2,1)),(Sys.sleep())) Note. Why I am using as.numeric is as I have a list of time at which I wanted to run the command. Something like 1349142243 1349138667 Question. 1. The query is running before the time reaches in both R-Studio and R-Terminal in Mac based System. 2. Sys.sleep () is ineffective in putting R on sleep till the time comes. Bhupendrasinh Thakre [[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] Ifelse Execution
Ifelse is a vector function and is absolutely inappropriate for that use. Use if instead. Also, read the help for Sys.sleep... you need to tell it how long you want to sleep. You should compute how long that is from now and sleep that long. --- 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. Bhupendrasinh Thakre vickytha...@gmail.com wrote: Hi Everyone, I am trying to run a time based query and need some of your help. Not much of data or packages. Just a simple one. Query I am trying to execute. ifelse ((as.numeric(as.POSIXct(2012-10-01 20:38:00))), (rnorm(1,2,1)),(Sys.sleep())) Note. Why I am using as.numeric is as I have a list of time at which I wanted to run the command. Something like 1349142243 1349138667 Question. 1. The query is running before the time reaches in both R-Studio and R-Terminal in Mac based System. 2. Sys.sleep () is ineffective in putting R on sleep till the time comes. Bhupendrasinh Thakre [[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] Transform pairwise observations into a table
On Oct 1, 2012, at 2:30 PM, arun wrote: HI AHJ, No problem. One more way in addition to reshape() (Rui's suggestion) to get the same result. library(reshape) as.matrix(cast(melt(dat1,id=c(ind1,ind2)),ind1~ind2,value=value)) # 1 2 3 4 #1 1.000 0.250 0.125 0.5 #2 0.250 1.000 0.125 0.5 #3 0.125 0.125 1.000 0.5 #4 0.500 0.500 0.500 1.0 A.K. That looks a tad ... well, ... complicated. So perhaps these base-only solutions with tapply might be more accessible: Some of them do border on the whimsical, I will admit: with (dat1, tapply(coef, list(ind1,ind2), I)) with (dat1, tapply(coef, list(ind1,ind2), c)) with (dat1, tapply(coef, list(ind1,ind2), ^, 1)) with (dat1, tapply(coef, list(ind1,ind2), +, 0)) It is a specific response to the request for a `table`-like function tha twouldallow the application of other functions. Cnage the `1` to `2` in the third instance and you get the tabulated squares. And do not forget the availability of `ftable` to flatten the output of `tapply` retunred values. -- David. - Original Message - From: Hadjixenofontos, Athena ahadjixenofon...@med.miami.edu To: arun smartpink...@yahoo.com Cc: R help r-help@r-project.org Sent: Monday, October 1, 2012 12:59 PM Subject: Re: [R] Transform pairwise observations into a table Thank you. I had looked at xtabs but misunderstood the syntax. This is great. :) AHJ On Oct 1, 2012, at 12:53 PM, arun smartpink...@yahoo.com wrote: Hi, Try this: dat1-read.table(text= ind1 ind2 coef 1 1 1 1 2 0.25 1 3 0.125 1 4 0.5 2 2 1 2 1 0.25 2 3 0.125 2 4 0.5 3 3 1 3 1 0.125 3 2 0.125 3 4 0.5 4 4 1 4 1 0.5 4 2 0.5 4 3 0.5 ,sep=,header=TRUE) mat1-as.matrix(xtabs(coef~ind1+ind2,data=dat1)) #ind2 #ind1 1 2 3 4 # 1 1.000 0.250 0.125 0.500 #2 0.250 1.000 0.125 0.500 #3 0.125 0.125 1.000 0.500 #4 0.500 0.500 0.500 1.000 A.K. - Original Message - From: AHJ ahadjixenofon...@med.miami.edu To: r-help@r-project.org Cc: Sent: Monday, October 1, 2012 12:17 PM Subject: [R] Transform pairwise observations into a table Hi, I have a table of pairs of individuals and a coefficient that belongs to the pair: ind1ind2coef 111 120.25 130.125 140.5 221 210.25 230.125 240.5 331 310.125 320.125 340.5 441 410.5 420.5 430.5 And I want to convert it to a matrix where each individual is both a row and a column and at the intersection of each pair is the coefficient that belongs to that pair: 1 2 3 4 11 0.25 0.1250.5 20.25 1 0.1250.5 30.1250.12510.5 40.5 0.5 0.51 If table() would allow me to specify something other than frequencies to fill the table with, it would be what I need. I tried a few different combinations of t() and unique() but none of it made enough sense to post as my starting code... I am just lost. Any help would be greatly appreciated. Thank you, AHJ -- View this message in context: http://r.789695.n4.nabble.com/Transform-pairwise-observations-into-a-table-tp4644706.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. David Winsemius, MD Alameda, CA, USA __ 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] nlme: spatial autocorrelation on a sphere
Please be assured, I really did RTFM (i.e. Pinheiro Bates) as well as all the online sources, and have accessed the codes of the corStruct functions, but they do not reveal how the metric argument is used. Therefore, I can't reprogram to include haversine distance, as seen in the ramps package. Smith et al. (2008) write The spatial correlation structures in nlme are not directly used because they do not allow great circle distance, which is very commonly needed for spatial data so there may be a case for adding this functionality to nlme. If anyone knows a way in to the relevant code, please let me know. Thanks Dan Reference Smith BJ, Yan J Cowles MK 2008. Unified Geostatistical Modeling for Data Fusion and Spatial Heteroskedasticity with R Package ramps. Journal of Statistical Software 25(10) From: David Winsemius [dwinsem...@comcast.net] Sent: 01 October 2012 17:32 To: Spencer Graves Cc: Dan Bebber; r-help@r-project.org Subject: Re: [R] nlme: spatial autocorrelation on a sphere On Oct 1, 2012, at 8:34 AM, Spencer Graves wrote: On 10/1/2012 12:38 AM, David Winsemius wrote: snip LMCTVTFY: http://cran.r-project.org/web/views/Spatial.html What's LMCTVIFY? Please accept my apologies for attempting to be cute and I also apologize to Mr Bebber for reading his request far too quickly. I was trying to use (C)RAN (T)ask (V)iews as a verb. LMGTFY for LMCTVIY produced only LMGTFY, at least for me. There's space for it on Wikipedia (RTFM) after LMGTFY (which I found using LMGTFY Wikipedia). I don't think it warrants being enshrined. I might also have written: require(sos); findFn(metric spherical latitude longitude) might produce. In this instance that approach found the ramps::corRSpher function, which has a haversine metric, much more quickly than did my subsequent efforts with RSeek, which I conducted after I realized that Bebber had already made a good faith effort at identifying resources. -- David Winsemius, MD Alameda, CA, USA __ 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] nlme: spatial autocorrelation on a sphere
Hi, Dan: On 10/1/2012 8:07 PM, Dan Bebber wrote: Please be assured, I really did RTFM (i.e. Pinheiro Bates) as well as all the online sources, and have accessed the codes of the corStruct functions, but they do not reveal how the metric argument is used. Therefore, I can't reprogram to include haversine distance, as seen in the ramps package. Smith et al. (2008) write The spatial correlation structures in nlme are not directly used because they do not allow great circle distance, which is very commonly needed for spatial data so there may be a case for adding this functionality to nlme. If anyone knows a way in to the relevant code, please let me know. Thanks for the clarification and for your persistence in this issue. Others clearly have been asking for an enhancement of this nature. I assume you've checked and confirmed that ramps::corRSpher, as suggested by David Winsemius, will NOT do what you want. If so, I suggest you subscribe to R-sig-mixed-models (per r-project.org - Mailing Lists), and post your question, citing Smith et al., to there. When you do, I suggest you include cc: Douglas Bates ba...@stat.wisc.edu. I'm sorry I can't help more, but I hope you will soon find a solution to this problem. Best Wishes, Spencer Thanks Dan Reference Smith BJ, Yan J Cowles MK 2008. Unified Geostatistical Modeling for Data Fusion and Spatial Heteroskedasticity with R Package ramps. Journal of Statistical Software 25(10) From: David Winsemius [dwinsem...@comcast.net] Sent: 01 October 2012 17:32 To: Spencer Graves Cc: Dan Bebber; r-help@r-project.org Subject: Re: [R] nlme: spatial autocorrelation on a sphere On Oct 1, 2012, at 8:34 AM, Spencer Graves wrote: On 10/1/2012 12:38 AM, David Winsemius wrote: snip LMCTVTFY: http://cran.r-project.org/web/views/Spatial.html What's LMCTVIFY? Please accept my apologies for attempting to be cute and I also apologize to Mr Bebber for reading his request far too quickly. I was trying to use (C)RAN (T)ask (V)iews as a verb. LMGTFY for LMCTVIY produced only LMGTFY, at least for me. There's space for it on Wikipedia (RTFM) after LMGTFY (which I found using LMGTFY Wikipedia). I don't think it warrants being enshrined. I might also have written: require(sos); findFn(metric spherical latitude longitude) might produce. In this instance that approach found the ramps::corRSpher function, which has a haversine metric, much more quickly than did my subsequent efforts with RSeek, which I conducted after I realized that Bebber had already made a good faith effort at identifying resources. -- David Winsemius, MD Alameda, CA, USA -- Spencer Graves, PE, PhD President and Chief Technology Officer Structure Inspection and Monitoring, Inc. 751 Emerson Ct. San José, CA 95126 ph: 408-655-4567 web: www.structuremonitoring.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] Basic question about: - and method start with dot.
Dear list, When I read some source code, I find lot of place used symbol - , e.g. lastTime - newTime; What is the meaning here? Also, I find some method with the name start with dot, e.g. .RowStandardizeCentered = function(x) { div = sqrt( rowSums(x^2) ); div[ div == 0 ] = 1; return( x/div ); } What is the special meaning for the method name start with a dot? Thank you very much in advance. __ 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] nlme: spatial autocorrelation on a sphere
Hi Spencer, thanks, I'll try R-sig-mixed-models. I tried using the corRSpher function from ramps in gls but received the following error: Error in recalc.corSpatial(object[[i]], conLin) : Unknown spatial correlation class Will post back here if I am successful. Dan From: Spencer Graves [spencer.gra...@prodsyse.com] Sent: 02 October 2012 04:29 To: Dan Bebber Cc: David Winsemius; r-help@r-project.org Subject: Re: [R] nlme: spatial autocorrelation on a sphere Hi, Dan: On 10/1/2012 8:07 PM, Dan Bebber wrote: Please be assured, I really did RTFM (i.e. Pinheiro Bates) as well as all the online sources, and have accessed the codes of the corStruct functions, but they do not reveal how the metric argument is used. Therefore, I can't reprogram to include haversine distance, as seen in the ramps package. Smith et al. (2008) write The spatial correlation structures in nlme are not directly used because they do not allow great circle distance, which is very commonly needed for spatial data so there may be a case for adding this functionality to nlme. If anyone knows a way in to the relevant code, please let me know. Thanks for the clarification and for your persistence in this issue. Others clearly have been asking for an enhancement of this nature. I assume you've checked and confirmed that ramps::corRSpher, as suggested by David Winsemius, will NOT do what you want. If so, I suggest you subscribe to R-sig-mixed-models (per r-project.org - Mailing Lists), and post your question, citing Smith et al., to there. When you do, I suggest you include cc: Douglas Bates ba...@stat.wisc.edu. I'm sorry I can't help more, but I hope you will soon find a solution to this problem. Best Wishes, Spencer Thanks Dan Reference Smith BJ, Yan J Cowles MK 2008. Unified Geostatistical Modeling for Data Fusion and Spatial Heteroskedasticity with R Package ramps. Journal of Statistical Software 25(10) From: David Winsemius [dwinsem...@comcast.net] Sent: 01 October 2012 17:32 To: Spencer Graves Cc: Dan Bebber; r-help@r-project.org Subject: Re: [R] nlme: spatial autocorrelation on a sphere On Oct 1, 2012, at 8:34 AM, Spencer Graves wrote: On 10/1/2012 12:38 AM, David Winsemius wrote: snip LMCTVTFY: http://cran.r-project.org/web/views/Spatial.html What's LMCTVIFY? Please accept my apologies for attempting to be cute and I also apologize to Mr Bebber for reading his request far too quickly. I was trying to use (C)RAN (T)ask (V)iews as a verb. LMGTFY for LMCTVIY produced only LMGTFY, at least for me. There's space for it on Wikipedia (RTFM) after LMGTFY (which I found using LMGTFY Wikipedia). I don't think it warrants being enshrined. I might also have written: require(sos); findFn(metric spherical latitude longitude) might produce. In this instance that approach found the ramps::corRSpher function, which has a haversine metric, much more quickly than did my subsequent efforts with RSeek, which I conducted after I realized that Bebber had already made a good faith effort at identifying resources. -- David Winsemius, MD Alameda, CA, USA __ 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.