Re: [R] German Map in package maps
On Wed, 23 May 2007, Maja Schröter wrote: > Hello everybody, > > I'm an absolute newbe in R, so please be gentle to me. > > I am looking for a german map in the package maps by Becker and Wilks. > > After that I want to draw the German cities from world.cities of this > package in the map and want to draw further a line from Berlin to Munich. > > Has anybody an idea? > There is no germany-only map in the maps databases, but the following might be close to what you want: > library(mapdata) # the world database in the maps package # probably isn't high enough resolution > map("worldHires", "Germany") > map.cities(country="Germany", minpop=1e6, capital=1) > lines(world.cities[world.cities$name %in% c("Berlin", "Munich") & world.cities$country=="Germany", c("long", "lat")]) > HTH Ray Brownrigg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] User defined grid on map
On Thu, 12 Apr 2007, Adriana Bejarano wrote: > This may be a trivial question for some of you: Is there a way to add > defined grid sizes on specific areas of a map?. I want to add 10kmx10km > grids to all California coastal counties: > > "california,alameda", "california,contra costa", "california,del norte", > "california,humboldt", "california,los angeles", "california,marin", > "california,mendocino" ,"california,monterey", "california,napa", > "california,orange", "california,riverside", "california,sacramento", > "california,san benito", "california,san bernardino", "california,san > diego" , "california,san francisco", "california,san joaquin", > "california,san luis obispo", "california,san mateo", "california,santa > barbara", "california,santa clara", "california,santa cruz", > "california,siskiyou", "california,solano", "california,sonoma", > "california,sutter", "california,trinity", "california,ventura", > "california,yolo" > The mapproj package has a map.grid() function, but that only allows you to specify latitude- and longitude-based grids. I don't think it would be too difficult to roll your own using standard trigonometry if you weren't too worried about the non-sphericity of the earth. Otherwise I guess there would be something 'out there' to provide you with some conversion between lat/long and distance from a fixed point. HTH, Ray Brownrigg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] How can i do the same thing in the China map?
usstata wrote: > The "maps" package has a function called "match.map", which is for map > coloring . > Its example is followed: > > # filled map showing Republican vote in 1900 > # (figure 6 in the reference) > data(state, package = "datasets") > data(votes.repub) > state.to.map <- match.map("state", state.name) > x <- votes.repub[state.to.map, "1900"] > gray.colors <- function(n) gray(rev(0:(n - 1))/n) > color <- gray.colors(100)[floor(x)] > map("state", fill = TRUE, col = color); map("state", add = TRUE) > > I want to do the same thing in the China map, but I can't find the Provinces > name of China. > Who can help me ? > nobody > > > a rookie > Tell us who you are, and you may get a more substantial reply. Ray Brownrigg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] circle fill problem
I suspect the number is: trunc((x + 3*r)/(2*r)) * trunc((y + 3*r)/(2*r)) + trunc((x + 2*r)/(2*r)) * trunc((y + 2*r)/(2*r)) where x and y are the dimensions of the rectangle and r is the radius of the circle. I have some code which graphs the overlapping circles if you want visually to check the algorithm (but although it is only a page, I won't post it here). HTH Ray Brownrigg On Thursday 08 February 2007 22:39, MINI GHOSH wrote: > Dear Ingmar and Robin, > > Thanks for you suggestions. I will see to it. > > Regards, > Mini > > --- Ingmar Visser <[EMAIL PROTECTED]> wrote: > > Robin & Mini, > > For those interested, googling for the 'orange > > packing problem' as it > > is known, or more officially the sphere packing > > problems gives you > > quite a few hits on these and similar problems. > > So at least the 3-d case the problem has been solved > > (I imagine the > > problem is easier in 2-d ...) > > hth, Ingmar > > > > On 8 Feb 2007, at 09:52, Robin Hankin wrote: > > > Mini > > > > > > This is a hard problem in general. > > > > > > Recreational mathematics has wrestled with > > > this and similar problems over the years; the > > > general field is the "set cover problem" but > > > in your case the sets are uncountably infinite > > > (and there are uncountably many of them). > > > > > > I would be surprised if your problem were not NP > > > > complete. > > > > > HTH > > > > > > > > > Robin > > > > > > On 8 Feb 2007, at 05:15, MINI GHOSH wrote: > > >> Dear R user, > > >> > > >> I want to know is there a way to find the minimum > > >> number of circles (of given radius) required to > > > > fill a > > > > >> given area (say rectangular) where overlapping of > > >> circles is allowed. > > >> > > >> Thanks, > > >> Regards, > > >> Mini Ghosh > > >> > > >> __ > > >> R-help@stat.math.ethz.ch mailing list > > >> https://stat.ethz.ch/mailman/listinfo/r-help > > >> PLEASE do read the posting guide > > > > http://www.R-project.org/posting- > > > > >> guide.html > > >> and provide commented, minimal, self-contained, > > > > reproducible code. > > > > > -- > > > Robin Hankin > > > Uncertainty Analyst > > > National Oceanography Centre, Southampton > > > European Way, Southampton SO14 3ZH, UK > > > tel 023-8059-7743 > > > > > > __ > > > R-help@stat.math.ethz.ch mailing list > > > https://stat.ethz.ch/mailman/listinfo/r-help > > > PLEASE do read the posting guide > > > > http://www.R-project.org/posting- > > > > > guide.html > > > and provide commented, minimal, self-contained, > > > > reproducible code. > > __ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide > http://www.R-project.org/posting-guide.html and provide commented, minimal, > self-contained, reproducible code. __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Cannot scale my map
On Tuesday 12 December 2006 09:07, Mihalis Tsoukalos wrote: > Dear list, > > I have the following problem: > > I want to plot my data and display a map at the same graph. The > problem is that my map is not properly sized (is very small). It > appears at the bottom-left side of the output. > I am using the following R commands: > > { > plot(Time, col="lightgrey") > title("Connections per Minute") > map("world", fill=TRUE, lwd=0.25, col="green", add=TRUE) > box() > } > > Please note that Time contains data in the following format: > 00:00 > 00:01 > ... > 22:04, etc. > > Can you please help me? > Well, the only way a map(..., add=TRUE) will work is if the original scale is in the order of latitude and longitude in degrees. It is not clear exactly how the data you are plotting relates to the map you want. If you just want the map as a 'background', then you will need to rescale one or other of the plots. Hope this helps, Ray Brownrigg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] indexing question
On Tuesday 14 November 2006 09:28, Leeds, Mark (IED) wrote: > thanks beilton but that won't work. A diff will also include 61-41 etc > and I don't want to include those. > > I'm working on using lapply or sapply with a seq along 31, 61, etc. > I'll let you know if it works. > Try looking at: dim(idx) <- c(3, length(idx)/3) then: tmp[idx[2, ]] - temp[idx[1, ]] temp[idx[3, ]] - temp[idx[2, ]] HTH Ray Brownrigg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] RE : avoiding a loop: "cumsum-like"
Well, I do have a solution which works for the data set you provide, but possibly not in a more general case. Firstly, tidying up your code, but using essentially the same looping algorithm, can provide a speed improvement of approximately 3:1. Here is a first attempt: mycode1 <- function(tab) { len <- diff(range(tab$Date)) + 1 res <- numeric(len) val <- 0 for (i in 1:len) { if (is.na(tab$posit.lat[i])) { val <- val + tab$x.jour[i] } else { if (res[tab$posit.lat[i]] < 30) { val <- val + tab$x.jour[i] } else { val <- val + tab$x.jour[i] + 0.8*res[tab$posit.lat[i]] } } res[i] <- val } return(res) } Then using a cumsum()-based algorithm can provide an overall 10:1 speed improvement: mycode2 <- function(tab) { res0 <- cumsum(tab$x.jour) res1 <- cumsum(ifelse(is.na(tab$posit.lat), 0, 0.8* (res0[tab$posit.lat] >= 30) * res0[tab$posit.lat])) res2 <- cumsum(ifelse(is.na(tab$posit.lat), 0, 0.8*res1[tab$posit.lat])) return(res0 + res1 + res2) } The condition is that: res[tab$posit.lat[tab$posit.lat[tab$posit.lat[length(tab$posit.lat) < 30 where tab is the data and res is the result. [There is also an implicit assumption that the result is monotonic.] HTH Ray Brownrigg GOUACHE David wrote: > Thanks Petr for taking a stab at it. > I have yet to figure out a way to do it, but if I do I'll post it. > Cheers > > David > > -Message d'origine- > De : Petr Pikal [mailto:[EMAIL PROTECTED] > Envoyé : vendredi 3 novembre 2006 09:05 > À : GOUACHE David; r-help@stat.math.ethz.ch > Objet : Re: [R] avoiding a loop: "cumsum-like" > > Hi > > I have not seen any answer yet so I wil try (partly). > > I believe that the loop can be vectorised but I am a little bit lost > in your fors and ifs. I found that first part of res is same as > cumsum(tab$x.jour) until about 81st value. However I did not decipher > how to compute the remaining part. I tried to add > cumsum(tab$posit.lat) (after changing NA to 0) what is not correct. > > Probably some combination of logical operation and summing can do > what you want. I thought that something like > ((cumsum(tab$posit.lat)*0.8)*(cumsum(tab$x.jour)>30)+cumsum(tab$x.jour > )) > > can do it but the result is defferent from your computation. > Not much of help, but maybe you can do better with above suggestion. > > Petr > > > > On 2 Nov 2006 at 11:15, GOUACHE David wrote: > > Date sent:Thu, 2 Nov 2006 11:15:49 +0100 > From: "GOUACHE David" <[EMAIL PROTECTED]> > To: > Subject: [R] avoiding a loop: "cumsum-like" > > >> Hello Rhelpers, >> >> I need to run the following loop over a large number of data-sets, and >> was wondering if it could somehow be vectorized. It's more or less a >> cumulative sum, but slightly more complex. Here's the code, and an >> example dataset (called tab in my code) follows. Thanks in advance for >> any suggestions! >> >> res<-0 >> for (i in min(tab$Date):max(tab$Date)) >> { >> if (is.na(tab$posit.lat[tab$Date==i])==T) >> { >> res<-c(res,res[length(res)]+tab$x.jour[tab$Date==i]) >> } >> else >> { >> if (res[tab$posit.lat[tab$Date==i]+1]<30) >> { >>res<-c(res,res[length(res)]+tab$x.jour[tab$Date==i]) >> } >> else >> { >>res<-c(res,res[length(res)]+tab$x.jour[tab$Date==i]+0.8*res[tab$pos >>it.lat[tab$Date==i]+1]) >> } >> } >> } >> res[-1] >> >> >> Date x.jour posit.lat >> 358040 NA >> 358050 NA >> 358060 NA >> 358070 NA >> 358080 NA >> 358092.97338883 NA >> 358102.796389915 NA >> 358110 NA >> 358120 NA >> 358131.000711886 NA >> 358140.894422571 NA >> 358150 NA >> 358160 NA >> 358170 NA >> 358180 NA >> 358190 NA >> 358200 NA >> 358210 NA >> 358220 NA >> 358230 NA >> 358240 NA >> 358250 NA >> 358260 NA >> 358270 NA >> 358280 NA >> 358290 NA >> 358300 NA >> 358310 NA >> 358320 NA >> 358330 NA >> 358340 NA >> 35
Re: [R] avoiding a loop
On Wednesday 25 October 2006 07:36, Leeds, Mark (IED) wrote: > I think I asked a similar question 3 years ago to the Splus list and I > think the answer was no or noone answered so noone should spend more > than 5 minutes on this > because it could definitely be a waste of time. > > My question is whether the function below can be rewritten without a for > loop. apply is fine if it can be done that way but i doubt it. I call it > a lot and would > prefer to not loop. > > #--- > -- > > constructLt<-function(invector) { > > outvector<-invector > > for ( i in 2:length(invector) ) { > if ( invector[i] < 1 ) { > outvector[i]<-invector[i]*outvector[i-1] > } > } > > return(outvector) > > } > Depending on the nature of your data, there is a faster way. It still involves looping, but not over the entire vector. Try the following: constructLt <- function(invector) { outvector <- invector cs <- cumsum(rle(invector < 1)$lengths) if (invector[1] < 1) cs <- c(1, cs) for (i in 0:(length(cs)%/%2 - 1)){ starti <- cs[2*i + 1] stopi <- cs[starti + 1] outvector[starti:stopi] <- cumprod(invector[starti:stopi]) } return(outvector) } It is in the order of 3 times as fast for random vectors of considerable length (> 1000). For random vectors of length 50 it is about the same speed as the full looping algorithm. However if the data is such that there are longer runs (than N(1, 1)), then you might expect a better speedup. HTH Ray Brownrigg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Possible bug?
On Thursday 05 October 2006 00:14, ronggui wrote: > Yes, the installer missed some file(s) in RHOME/bin. > When I choose "Message translations", the /bin directory contains the > following files. > -- snip -- I meant to mention that when I used the installer I didn't have to choose "Message translations", it was ticked by default. HTH, Ray Brownrigg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Possible bug?
On Wednesday 04 October 2006 14:17, Duncan Murdoch wrote: > On 10/3/2006 9:00 PM, ronggui wrote: > > This morning I downloaded R-2.4.0 and install in under Windows. I > > customized the installation and choose "Message translations",but I could > > not launch Rgui.exe successfully( Rterm.exe worked fine). If I did't > > choose "Message translations", Rgui.exe worked fine. > > > > "Message translations" is Simplified Chinese. > > Did you test any of the betas or release candidates? I run an English > language version of Windows, and I don't even get offered the chance to > install in Simplified Chinese, so I certainly didn't test this. > > Duncan Murdoch > > 7*runif(1) [1] 2.897160 Well, I have a (self compiled) R-2.4.0 beta (2006-09-24 r39497) which doesn't exhibit the behaviour in the same environment (Windows XP Professional English version but with Simplified Chinese set as the locale) that the released R-2.4.0-win32.exe does. [Perhaps my version didn't set up the equivalent of "Message translations" so this may be a red herring, but it does output the startup message in Chinese, and error messages are also output in Chinese.] What happens with the release version is an error message: "R for Windows GUI front-end has encountered a problem and needs to close. We are sorry for the inconvenience." Then it offers to "Please tell Microsoft about the problem", and if you "click here" the error signature has: AppName: rgui.exe AppVer: 2.40.39566.0ModName: r.dll ModVer: 2.40.39566.0Offset: 000f22b3 There is a lot more information as well, but I don't know how much is relevant. By "Simplified Chinese set as the locale" I mean 'Control Panel/Regional and Language Options' with Chinese (PRC) set in 'Regional Options/Standards and Formats' and 'Advanced/Language for non-Unicode Programs'. [I don't rely on the Chinese environment, so this doesn't materially affect me, but I may be able to help with diagnosis.] Hope this helps, Ray Brownrigg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Speeding indexing and sub-sectioning of 3d array
> Date: Thu, 10 Aug 2006 14:34:27 -0400 > From: "Swidan, Firas" <[EMAIL PROTECTED]> > > Hi Patrick, > > Thanks for the help. The function I listed is just an example. I isolated > and kept only the problematic part in my code for clarity sake. I ended up > implementing the functionality in C and now it takes 22 seconds to calculate > the objective. > > Best regards, > Firas. > Interestingly, I was able to develop an algorithm in R that achieves the same order-of-magnitude speedup as your C code, but at the expense of greater memory requirements. However it only works if the function you are using is really is mean() [your code labels use Median]. It does this by making use of cumsum() and logical indexing, working with sums of values rather than calculationg the means and then dividing by the numbers of values in the hypercube at the end. If you want to try coding this algorithm in C for even greater performance improvement (or for interest only), let me know. I suspect it will be difficult to code in C because of the vectorisation it takes advantage of. In the output below, cK3d() is your algorithm (slightly adjusted to cover the whole matrix and to return something), and cK3dme is my equivalent, running on a Pentium IV 3.2GHz, NetBSD system with 1GB memory. Regards, Ray Brownrigg > x <- rnorm(245*175*150) > dim(x) <- c(245, 175, 150) > unix.time(yme <- cK3dme(x, 3)) [1] 13.870 1.690 15.813 0.000 0.000 > unix.time(y <- cK3d(x, 3)) [1] 500.206 0.035 505.738 0.000 0.000 > all.equal(y, yme) [1] TRUE > __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] add city and point in the map
> From: YIHSU CHEN <[EMAIL PROTECTED]> > Date: Thu, 27 Apr 2006 10:50:49 -0400 > > I'm trying to use packages "maps" and "mapdata" (see blow) to display > the research resutls on map (Mid-Atlantic region). In particular, I > need to mark a number of points in the map by giving their latitude and > longitude information. For instance, I would like to mark a point on > (long, lat) =(75.56027, 39.09271). Also, I need to mark several cities > that I suspect might not in the "world.cities" database such as > Baltimore, etc. Presumably, I can do the same thing by finding out the > lat and long infor of Baltimroe and add text next to the point I mark. > > Thanks a billion > > Yihsu > > > library(maps) > library(mapdata) > state.list <- c("maryland","new jersey","west virginia", > "virginia","ohio","delaware","pennsylvania","kentucky","indiana","north > carlina") > map("state",state.list) > map.cities(x=world.cities, Country="US",minpop=10) > Several things: 1) the US is west of longitude 0, so you need longitude -75.56027 2) map.cities uses x=world.cities by default, so you don't need to specify it in your call, but if you do specify it, you need to precede it with a call to: data(world.cities) [I guess this is a bug, perhaps just at the documentation level] 3) map.cities takes country= as an argument, not Country= 4) The US is identified in the maps package as USA, so you would use country="USA" (but you don't need country= anyway, beacause your map is only of the USA, so there are no cities to exclude). 5) Baltimore is in world.cities, try map.cities(minpop=50, label=T) 6) If you really want North Carolina, you need to spell it correctly So what you need is: library(maps) # library(mapdata) # not needed for this state.list <- c("maryland", "new jersey", "west virginia", "virginia", "ohio", "delaware", "pennsylvania", "kentucky", "indiana", "north carolina") map("state", state.list) map.cities(minpop=10) Hope this helps, Ray Brownrigg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] matching a given sd range
> Date: Fri, 31 Mar 2006 14:44:05 -0800 (PST) > From: "Fred J." <[EMAIL PROTECTED]> > > given a numeric array (a sequence of reals), I am interested in > finding the subsets of sequences (each with start and end index) which match > a given sd range. > > I read the docs on match and which and the "see also" but could not come up > with a way. I could loop with a stepping window over the sequence but that > would be limited to a fixed size window, I guess I as well loop through all > window sizes from 1:length(data) but that is not elegant. any hints towards > a solution? > Presumably your "array" is a vector, in which case you can use vectorisation to eliminate one of the loops. The following function produces a matrix in upper triangular form for which res[i, j] is the variance of the sequence dat[i:j] (where dat is your "array"). It is (approximately) two orders of magnitude faster than the double loop, but doesn't produce exactly the same numbers because of numerical rounding because it uses the formula: variance = (sum(xi^2) - n*xbar^2)/(n - 1) Now all you have to do is test this matrix against the range of variances you want to isolate, and use the matching indices to provide you with the endpoints of your subsets. Hope this helps, Ray Brownrigg allvar <- function(dat) { len <- length(dat) res <- numeric((len - 1)*len/2) mym <- dat[1] myv <- 0 for (j in 2:len) { oldm <- mym mym <- (((j - 1):1)*mym + dat[j])/(j:2) myv <- (((j - 2):0)*myv + ((j - 1):1)*oldm^2 + dat[j]^2 - (j:2)*mym^2)/ ((j - 1):1) res[((j - 2)*(j - 1)/2 + 1):((j - 1)*j/2)] <- myv mym <-c(mym, dat[j]) myv <- c(myv, 0) } return(res) } # Now test it > n <- 500; dat <- runif(n) > unix.time({res2 <- matrix(0, n, n); res2[upper.tri(res2)] <- allvar(dat)}) [1] 0.64 0.00 0.64 0.00 0.00 > unix.time({res1 <- matrix(0, n, n); for (j in 2:n) for (i in 1:(j - 1)) > res1[i, j] <- var(dat[i:j])}) [1] 72.64 13.74 86.38 0.00 0.00 > max(abs(res1 - res2)) [1] 1.249001e-15 > __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Plotting shapefiles on existing maps
> From: [EMAIL PROTECTED] > > I have plotted a map of the Barents Sea and surrounding coastline using: > > map('worldHires',ylim=c(50,85),xlim=c(5,65),fill=T,resolution=0) > map.axes() > map.scale(x=30,metric=T) > > Next, I imported a shapefile with depth contours for the sea: > > contours<-read.shape("D://My Documents/BarentsSea.shp",dbf.data=T) > > (This is in mercator projection). > > Despite extensive searches of the help files and R site, I cannot find a > way to plot the contours onto the map. Does anyone have any suggestions? > It is not clear to me exactly what you are having difficulty with, but if you are able to produce a contour plot from that data, then note that both contour() and map() have an "add" argument, so as long as the axis limits are conformant (i.e. both longitude/latitude) you can add one plot to the other. Perhaps the trick is in converting your contours data into longitude/latitude form for the x and y components. Hope this helps, Ray Brownrigg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] bug in map('world') ?
> From: Joerg van den Hoff <[EMAIL PROTECTED]> > > hi, > > did'nt see anything in the archive: > > map('world',pro='rectangular',para=0) > > yields a strange artifact (horizontal bar) extending over the whole map > at a certain latitude range (approx 65 deg. north), whereas > > map('world',pro='rectangular',para=180) > > (which should be the same) does not show the artifact. > > the artifact shows up in other projections as well, e.g. > > map('world',pro='bonne',para=45) > > > which seems to show that the problem might be in the data base of the > map (i.e. polygon definition)?? > > > any ideas? > The short answer is RTFM, and now that I have again :-), I note that there is a wrap= option. Setting that to TRUE seems to have the desired effect. Ray __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] bug in map('world') ?
I seem to recall I came across something similar recently (well, relatively recently, perhaps a couple of years ago). The problem is related to how the code handles polygons that are split across the map boundaries. As I recall, the fix was to modify the polygons that straddle the date line. I'll have to delve into this again, but don't hold your breath, sorry. Ray Brownrigg > From: Joerg van den Hoff <[EMAIL PROTECTED]> > > hi, > > did'nt see anything in the archive: > > map('world',pro='rectangular',para=0) > > yields a strange artifact (horizontal bar) extending over the whole map > at a certain latitude range (approx 65 deg. north), whereas > > map('world',pro='rectangular',para=180) > > (which should be the same) does not show the artifact. > > the artifact shows up in other projections as well, e.g. > > map('world',pro='bonne',para=45) > > > which seems to show that the problem might be in the data base of the > map (i.e. polygon definition)?? > > > any ideas? > > regards, > > joerg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Installing packages without clicking
> I need to install several (too many) packages from local *.zip > files. > Are you aware of being able to select multiple .zip files at one time using SHIFT-click or CTRL-click? > is there any form to do it without clicking? > > I'm asking for a R code that allow me perform these task. > You could just do what the clicking does, i.e. call the function install.packages() with appropriate parameters. See utils:::menuInstallLocal Hope this helps, Ray Brownrigg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] expand.grid without expanding
> From: =?iso-8859-1?q?Lu=EDs_Torgo?= <[EMAIL PROTECTED]> > Date: Wed, 8 Feb 2006 18:08:40 + > > Dear list, > I've recently came across a problem that I think I've solved and that I > wanted > to share with you for two reasons: > - Maybe others come across the same problem. > - Maybe someone has a much simpler solution that wants to share with me ;-) > > The problem is as follows: expand.grid() allows you to generate a data.frame > with all combinations of a set of values, e.g.: > > expand.grid(par1=-1:1,par2=c('a','b')) > par1 par2 > 1 -1a > 20a > 31a > 4 -1b > 50b > 61b > > There is nothing wrong with this nice function except when you have too many > combinations to fit in your computer memory, and that was my problem: I > wanted to do something for each combination of a set of variants, but this > set was to large for storing in memory in a data.frame generated by > expand.grid. A possible solution would be to have a set of nested for() > cycles but I preferred a solution that involved a single for() cycle going > from 1 to the number of combinations and then at each iteration having some > form of generating the combination "i". And this was the "real problem": how > to generate a function that picks the same style of arguments as > expand.grid() and provides me with the values corresponding to line "i" of > the data frame that would have been created bu expand.grid(). For instance, > if I wanted the line 4 of the above call to expand.grid() I should get the > same as doing: > > expand.grid(par1=-1:1,par2=c('a','b'))[4,] > par1 par2 > 4 -1b > > but obviously without having to use expand.grid() as that involves generating > a data frame that in my case wouldn't fit in the memory of my computer. > > Now, the function I've created was the following: > > getVariant <- function(id,vars) { > if (!is.list(vars)) stop('vars needs to be a list!') > nv <- length(vars) > lims <- sapply(vars,length) > if (id > prod(lims)) stop('id above the number of combinations!') > res <- vector("list",nv) > for(i in nv:2) { > f <- prod(lims[1:(i-1)]) > res[[i]] <- vars[[i]][ceiling(id / f)] > id <- id - (ceiling(id/f)-1)*f > } > res[[1]] <- vars[[1]][id] > names(res) <- names(vars) > res > } > -- > > expand.grid(par1=-1:1,par2=c('a','b'))[4,] > par1 par2 > 4 -1b > > getVariant(4,list(par1=-1:1,par2=c('a','b'))) > $par1 > [1] -1 > > $par2 > [1] "b" > > I would be glad to know if somebody came across the same problem and has a > better suggestion on how to solve this. > A few minor improvements: 1) let id be a vector of indices 2) use %% and %/% instead of ceiling (perhaps debateable) 3) return a data frame as does expand.grid So your function now looks like: getVariant <- function(id, vars) { if (!is.list(vars)) stop('vars needs to be a list!') nv <- length(vars) lims <- sapply(vars, length) if (any(id > prod(lims))) stop('id above the number of combinations!') res <- vector("list", nv) for(i in nv:2) { f <- prod(lims[1:(i-1)]) res[[i]] <- vars[[i]][(id - 1)%/%f + 1] id <- (id - 1)%%f + 1 } res[[1]] <- vars[[1]][id] names(res) <- names(vars) return(as.data.frame(res)) } Now, for example, you get: > expand.grid(par1=-1:1,par2=c('a','b'),par3=c('w','x','y','z'))[12:15,] par1 par2 par3 121bx 13 -1ay 140ay 151ay > getVariant(12:15,list(par1=-1:1,par2=c('a','b'), par3=c('w','x','y','z'))) par1 par2 par3 11bx 2 -1ay 30ay 41ay > Note that you will run into trouble when the product of the lengths is greater than the largest representable integer on your system. Hope this helps, Ray Brownrigg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Number of replications of a term
There's an even faster one, which nobody seems to have mentioned yet: rep(l <- rle(ids)$lengths, l) Timing on my 2.8GHz NetBSD system shows: > length(ids) [1] 45150 > # Gabor: > system.time(for (i in 1:100) ave(as.numeric(factor(ids)), ids, FUN = length)) [1] 3.45 0.06 3.54 0.00 0.00 > # Barry (and others I think): > system.time(for (i in 1:100) table(ids)[ids]) [1] 2.13 0.05 2.20 0.00 0.00 > Me: > system.time(for (i in 1:100) rep(l <- rle(ids)$lengths, l)) [1] 1.60 0.00 1.62 0.00 0.00 Of course the difference between 21 milliseconds and 16 milliseconds is not great, unless you are doing this a lot. Ray Brownrigg > From: Gabor Grothendieck <[EMAIL PROTECTED]> > > Nice. I timed it and its much faster than mine too. > > On 1/24/06, Barry Rowlingson <[EMAIL PROTECTED]> wrote: > > Laetitia Marisa wrote: > > > Hello, > > > > > > Is there a simple and fast function that returns a vector of the number > > > of replications for each object of a vector ? > > > For example : > > > I have a vector of IDs : > > > ids <- c( "ID1", "ID2", "ID2", "ID3", "ID3","ID3", "ID5") > > > > > > I want the function returns the following vector where each term is the > > > number of replicates for the given id : > > > c( 1, 2, 2, 3,3,3,1 ) > > > > One-liner: > > > > > table(ids)[ids] > > ids > > ID1 ID2 ID2 ID3 ID3 ID3 ID5 > > 1 2 2 3 3 3 1 > > > > 'table(ids)' computes the counts, then the subscripting [ids] looks it > > all up. > > > > Now try it on your 40,000-long vector! > > > > Barry __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] elements in a matrix to a vector
> my apologies to everyone taking the time to answer this question. i > didn't explain myself very well. > > i have the indices of the elements of the matrix and the matrix won't > necessarily have zeros everywhere else. is there a way to use the > indices to strip the elements i need? > RTFM ?"[". > x <- matrix(c(0, 2, 0, 0, 0, 4, 3, 0, 0), nrow=3) > x [,1] [,2] [,3] [1,]003 [2,]200 [3,]040 > ?"[" > xx <- matrix(c(1, 2, 3, 3, 1, 2), ncol=2) > xx [,1] [,2] [1,]13 [2,]21 [3,]32 > x[xx] [1] 3 2 4 > HTH, Ray __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Error in integrate
Murray Jorgensen <[EMAIL PROTECTED]> wrote: > I'm using R 2.0.1 under Windows XP. I get the following message when I > run the code listed below. I don't seem to have any problems using the > function "slice" outside "integrate". > > Error in integrate(slice, 0, Inf, , , , , , , a, b) * delta : > non-numeric argument to binary operator > RTFM! Integrate returns: Value: A list of class '"integrate"' with components : : You need: A[j] <- integrate(slice,0,Inf,,,a,b)$value*delta Ray __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] generating a vector from clusters of logicals
> From: Troels Ring <[EMAIL PROTECTED]> > > I have a vector of clusters of TRUE and FALSE like > c(TRUE,TRUE,TRUE...,FALSE,FALSE,FALSE,TRUE,TRUE...) and want to make > that into a vector > of c(1,1,1,1...2,2,2,2,.3,3,3,3) increasing the number assigned to each > cluster as they change. > How would I do that ? > How about: > TF <- c(TRUE,TRUE,TRUE,FALSE,FALSE,FALSE,TRUE,TRUE,FALSE) > rep(1:length(rlel <- rle(TF)$lengths), rlel) [1] 1 1 1 2 2 2 3 3 4 HTH, Ray Brownrigg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
[R] problem with Hershey fonts
This was reported to me by a colleague in China, so I may not be reproducing exactly what they are seeing (which I suspect is rw2011), but this is what I see: > version _ platform i386--netbsdelf arch i386 os netbsdelf system i386, netbsdelf status major2 minor1.1 year 2005 month06 day 20 language R > help(Hershey) : : If the 'vfont' argument to one of the text-drawing functions ('text', 'mtext', 'title', 'axis', and 'contour') is a character vector of length 2, Hershey vector fonts are used to render the text. : The other useful escape sequences all begin with '\\'. These are described below. Remember that backslashes have to be doubled in R character strings, so they need to be entered with _four_ backslashes. : : > plot(runif(100)) # to get something on the screen > text(0, 1, '\\Re', vfont=c('serif', 'plain')) # works > title(main = '\\Re', vfont=c('serif','plain'))# doesn't work, and... Warning message: parameter "vfont" could not be set in high-level plot() function > mtext('\\Re', 2, vfont=c('serif','plain'))# doesn't work, but... Warning message: Hershey fonts not yet implemented for mtext() in: mtext(text, side, line, outer, at, adj, padj, cex, col, font, > axis(3, at=50, tick=F, labels='\\Re', vfont=c('serif','plain')) # > doesn't work Now at least mtext() 'redeems' itself with its Warning, but title() contradicts the help(Hershey), and axis() just gets it wrong by printing the characters "\Re" (or "\\Re" when the string is specified as 'Re' as alluded to in the Hershey documentation - but somewhat contradicted by the success of the text() call above). Now perhaps the problem is just with the Hershey documentation (or my reading of it), but I don't know what is supposed to work. Or is there some combination that I haven't tried? Ray Brownrigg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Maps, Eastern Europe
> From: "Rau, Roland" <[EMAIL PROTECTED]> > > I would like to employ a European map in a presentation. My idea was to > use: > library(mapdata) > map("worldHires", c("Austria", "Switzerland", "Germany")) > > where I include all countries from my analysis as a vector of character > strings like in the example above. Unfortunately, I was unable to > specify the Baltic States (Latvia, Lithuania, Estonia) or the Czech > Republic. > Initially, I thought I misspelled the country names. To avoid this > obvious error, I checked within 'mapdata_2.0-14.tar.gz' the file > 'worldHires.names'. To my surprise, the spelling was not the source of > the error - it was rather the case that those countries were not > included. > > Did some people experience similar problems? And if yes: how did you > proceed? > Unfortunately the world databases for the maps package are very old (early 1990's) and have not been updated since, except for minor typographical errors in country names. Note that Czechoslovakia is there (as is USSR!). Also, it currently satisfies my needs, which is why I don't have the resources to update it. I suggest you start with the maptools package, but as I understand it, you need to also find some appropriate shapefile data. [Other respondents may improve upon this suggestion.] Regards, Ray Brownrigg maps package maintainer __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] a question about box counting
I said: > myfun <- function(x, y, ints) { > fx <- x %/% (1/ints) > fy <- y %/% (1/ints) > txy <- hist(fx + ints*fy+ 1, breaks=0:(ints*ints), plot=FALSE)$counts > dim(fxy) <- c(ints, ints) ^^^ > return(txy) > } Of course it should be: dim(txy) <- c(ints, ints) ^^^ Sorry about that, Ray __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] a question about box counting
> From: Deepayan Sarkar <[EMAIL PROTECTED]> Mon, 4 Apr 2005 13:52:48 -0500 > > On Monday 04 April 2005 13:22, Rajarshi Guha wrote: > > Hi, > > I have a set of x,y data points and each data point lies between > > (0,0) and (1,1). Of this set I have selected all those that lie in > > the lower triangle (of the plot of these points). > > > > What I would like to do is to divide the region (0,0) to (1,1) into > > cells of say, side = 0.01 and then count the number of cells that > > contain a point. > > > > My first approach is to generate the coordinates of these cells and > > then loop over the point list to see whether a point lies in a cell > > or not. > > > > However this seems to be very inefficient esepcially since I will > > have 1000's of points. > > > > Has anybody dealt with this type of problem and are there routines to > > handle it? > > A combination of cut and table/xtabs should do it, e.g.: > > > x <- runif(3000) > y <- runif(3000) > > fx <- cut(x, breaks = seq(0, 1, length = 101)) > fy <- cut(y, breaks = seq(0, 1, length = 101)) > > txy <- xtabs(~ fx + fy) > : Another significantly faster way (but not generating row/column names) is: x <- runif(3000) y <- runif(3000) ints <- 100 myfun <- function(x, y, ints) { fx <- x %/% (1/ints) fy <- y %/% (1/ints) txy <- hist(fx + ints*fy+ 1, breaks=0:(ints*ints), plot=FALSE)$counts dim(fxy) <- c(ints, ints) return(txy) } myfun(x, y, ints) Hope this helps, Ray Brownrigg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Recursive default argument reference
> From: boston knot <[EMAIL PROTECTED]> > Date: Mon, 24 Jan 2005 20:50:53 -0800 (PST) > "Recursive default argument reference" keeps appearing when I try to run a > haplo.score function in R for Windows. I'm new to using this program. Does > anyone know what this means? > I don't know exactly what it means, but I have seen it quite some time ago when I was creating a package using delay(). It is easily reproduced by: > x <- delay({runtimeerror; 10}) > x Error: Object "runtimeerror" not found > x Error: recursive default argument reference > HTH, Ray Brownrigg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] zero index and lazy evaluation in ifelse()
> Date: Thu, 13 Jan 2005 12:42:59 -0800 > From: [EMAIL PROTECTED] > > I don't understand this behavior: > > > a <- c(0, 1, 2, 3) > > b <- c(1, 2, 3, 4) > > ifelse (a == 0, 0, b[a]) > [1] 0 2 3 1 > > rather than the desired 0 1 2 3. Thanks for any explanation. > Look at b[a]. Ray __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] load object
> From: Douglas Bates <[EMAIL PROTECTED]> Fri Jan 14 08:35:33 2005 > > Weiwei Shi wrote: > > Hi, > > I happen to re-write my codes to save memory and my > > approach is write my obj into file first and later I > > load it. > > > > However, it seems like: > > load(filename) can load the object but the function > > returns the name of the object instead of the > > reference to it. For example, I have an object called > > r0.prune, which is saved by > > save(r0.prune, file='r0.prune') > > > > and later, I want to load it by using: > > load('r0.prune') > > but I need to put the reference to the object r0.prune > > into a var or a list. I tried: > > t<-load('r0.prune'), > > and class(t) gave me a char, which means t stores the > > name of obj instead of the obj itself. > > > > Sorry for the dumb question but please help... > > > > Weiwei > > I was going to suggest that you read the help page for load but when I > looked at it myself I found that it was not obvious what the effect of > calling load at the prompt is. The help page is accurate but it is a > bit confusing if you don't know what the default environment is. > > Anyway, when called from the prompt, load has the side effect of loading > the object into the global environment. Try > Well ?load does say: Value: A character vector of the names of objects created, invisibly. Note the plurals. The point is that the file being loaded may contain the definition of more than one R object. You can say: tt <- get(load('r0.prune')) # t is not a good example name to use which will do what you want, but: 1) load() has a side-effect of also creating r0.prune (in your case) 2) if the file 'r0.prune' contains more than one object, only the first is assigned to tt (but all of them are loaded into memory with their original names). Hope this helps, Ray Brownrigg __ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] counting numbers without replicates in a vector
> I am just wondering if there is an easy way to count > in a numeric vector how many numbers don't have > replicates. > For example, > a=c(1,1,2,2,3,4,5), how can I know there are three > numbers (3, 4 and 5) without replicates? > How about: length(table(a)[table(a) == 1]) ? Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Strange error from R CMD INSTALL
> Date: Wed, 08 Dec 2004 09:33:22 -0400 > From: Kjetil Brinchmann Halvorsen <[EMAIL PROTECTED]> > > I am trying to install a local package and get this unexpected > error: > > -- Making package UMSA > adding build stamp to DESCRIPTION > installing R files > installing data files > installing man source files > installing indices > Error: couldn't find function "na.omit" > Execution halted > > na.omit of course is in package stats, and that is listed in the > Depends field in the DESCRIPTION file. > Check what is in your data directory. Does something in there use (implicitly) na.omit? I have seen something like this (since 2.0.0) where a .R file in the data/ directory used a function within the package, but that function was not 'available' at the "installing indices" phase of the check/build. Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Manipulating contour plots
From: Caroline Houldcroft <[EMAIL PROTECTED]> Date: 01 Dec 2004 15:16:55 + > I would like to create a filled.contour plot with colour bar and overlay > a map. > > I need to create several of these representing variables with different > spatial extents for comparison. When I try to fix the plot area size and > axes limits as follows: > > par(fig=c(0,1.0,0.3125,1.0) > > filled.contour(x=seq((320/2)-180.5,(480/2)-180.5,len=nrow(oct94)), > y=seq((190/2)-90.5,(240/2)-90.5,len=ncol(oct94)), > z=as.matrix(oct94), > color=heat.colors, > levels=seq(0,30,2), > plot.axes={axis(1,seq(-60,60,by=20)); > axis(2,seq(5,30,by=5))}) > > > my definitions were ignored producing a stretched image of the area with > limits defined by the contour plot. > I cannot answer this part of the question because your example is not reproducible (Error in nrow(oct94) : Object "oct94" not found) [and the first line contains a syntax error.] > > I have also had difficulties overlaying a map on the same axes. The > following example, taken from the r-help archives, gave the error > "couldn't find function "world"" > > library(maps) > filled.contour(x=seq((320/2)-180.5,(480/2)-180.5,len=nrow(oct94)), > y=seq((190/2)-90.5,(240/2)-90.5,len=ncol(oct94)), > z=as.matrix(oct94), > color=heat.colors, > levels=seq(0,30,2), > plot.axes={world(ylim=c(5,30),xlim=c(-30,50),add=T); > axis(1);axis(2)}) > The world() function is not part of the maps package, it is in the fields package. If you want to use the maps package, try: map("world", ylim=c(5,30), xlim=c(-30, 50), add=T) Hope this helps, Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] library(fields) world shift function not working anymore
> I just upgraded to 2.01 on Mac OS 10.3.6. I used to use the command > (on R 1.9.x): > > world(ylim=c(-30,30), xlim = c(0,360), shift=TRUE, add=TRUE) > > to draw a world outline over my image plots. My data uses longitude > from (0, 360) so I need to use the shift function. After I upgraded, I > get the following error: > > > world(ylim=c(-30,30), xlim = c(0,360), shift=TRUE, add=TRUE) > Error in world(ylim = c(-30, 30), xlim = c(0, 360), shift = TRUE, add = > TRUE) : > NAs are not allowed in subscripted assignments > > Does anyone know a workaround for this? > Well, a workaround would be: library(maps) map("world2", ylim=c(-30,30), xlim = c(0,360), add = TRUE) (at least until fields is updated). CHANGES IN R VERSION 2.0.0 o Subassignments involving NAs and with a replacement value of length > 1 are now disallowed. (They were handled inconsistently in R < 2.0.0, see PR#7210.) For data frames they are disallowed altogether, even for logical matrix indices (the only case which used to work). Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] IFELSE across large array?
> From: "Liaw, Andy" <[EMAIL PROTECTED]> > Date: Tue, 23 Nov 2004 12:28:48 -0500 > > I'll give it half a crack: > > Steps a through c can be done via nested ifelse(), treating A and M as > vectors (as they really are). Step d is the hard one. I'd write a simple > Fortran code and use .Fortran() for that. > > I don't see how any of the *apply() functions can help here, as your > operations are element-wise, not dimension-wise. > > Andy > The original message mentions that the value 10 is actually "NODATA", and if one recodes 10 as NA, steps a) to c) become trivial, namely: A[A == 10] <- NA M[M == 10] <- NA return(A + M - 1) Then if step d) is performed first (i.e. appropriate values in A are replaced by the 'most common neighbour' [perhaps using round(mean(.., na.rm=T))] this still works, but would have to be repeated for each replication (the third dimension). Ray Brownrigg > > From: Sander Oom > > > > Dear all, > > > > As our previous email did not get any response, we try again with a > > reformulated question! > > > > We are trying to do something which needs an efficient loop > > over a huge > > array, possibly functions such as apply and related (tapply, > > lapply...?), but can't really understand syntax and examples in > > practice...i.e. cant' make it work. > > > > to be more specific: > > we are trying to apply a mask to a 3D array. > > By this I mean that when "overlaying" [i.e. comparing element > > by element] > > the mask on to the array the mask should change array > > elements according to > > the values of both array and mask elements > > > > the mask has two values: 1 and 10. > > > > the array elements have 3 values: 0, 1, or 10 > > > > sticking for the moment to the single 2d array case > > > > for example: > > [A= array] 100 10 1 10 0 > > 1 10 1 0 0 10 > > > > [ M=mask] 1 10 10 1 1 1 > > 101 1 1 10 10 > > > > I would like the array elements to: > > > > a) IF A(ij) !=10 and Mij = 1 > > leave A(ij) unchanged > > > > b) IF A(ij) != 10 and M(ij) =10 > > change A(ij) to M(ij) i.e mask value (10) > > > > c)IF A(ij) = 10 and M(ij) = 10 > >leave (Aij) unchanged > > > > d) IF A(ij) = 10 and M(ij) !=10 > > replace A(ij) with the majority value in the > > 8-neighborhood > > > > (or whatever if it is an edge element) BUT ignoring 10s in this > > neighborhood (i.e. with either 1 or 0, whichever is in majority) > > > > because the array is 3d I would like to repeat the thing with > > all the k > > elements (2d sub-arrays) of the array in question, using the > > same mask for > > al k elements > > > > Would you be able to suggest a strategy to do this? > > > > thanks very much > > > > Alessandro and Sander. __ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Re: 3d Map with bars
> From: [EMAIL PROTECTED] > Date: Fri, 19 Nov 2004 16:53:12 -0500 > > Thanks for reply. I need to first draw the map of USA a perspective plot. > I guess thats where my problem was. > Try something like this: library(maps) states <- map("state", plot=F) x1 <- rep(0, 3) x2 <- rep(0, 3) maxz <- 1 z <- matrix(c(0, 0, 0, 0, 0, 0, maxz, maxz, maxz), 3, 3) x1[-2] <- states$range[1:2] x2[-2] <- states$range[3:4] x1[2] <- x1[3] - 1e-6 x2[2] <- x2[3] - 1e-6 pmat <- persp(x1, x2, z, xlab = "", ylab = "", zlab="", axes = F, theta=0, phi=20, d=10) lines(trans3d(states$x, states$y, 0, pmat)) latitude <- -90 longitude <- 37 myz <- 0.6 lines(trans3d(rep(latitude, 2), rep(longitude, 2), c(0, myz), pmat), col=2) where trans3d is as defined in the examples for persp(). Hope this helps, Ray Brownrigg > Uwe Ligges <[EMAIL PROTECTED]> > 11/19/2004 04:33 PM > > To: [EMAIL PROTECTED] > cc: [EMAIL PROTECTED] > Subject:Re: 3d Map with bars > > > [EMAIL PROTECTED] wrote: > > > Apologies in advance for the question. I am trying to draw a map of the > US > > as a surface plot so that I would be able to drop bars on the different > > states (something like Uwe Ligges' scatterplot3d example 4). I am not > sure > > where to start looking for such a beast. If anyone has any pointers, > > ideas, I will be grateful. > > > > TIA, > > Partha > > How to "drop bars" with persp() has been described on R-help yesterday > or today, please check the mailing list's archives. > __ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] adjusting the map of France to 1830
> Date: Fri, 19 Nov 2004 15:59:25 -0500 > From: Michael Friendly <[EMAIL PROTECTED]> > > Here's what I tried. I can plot a selection of regions, but I > can't seem to remove an arbitrary list of region numbers, unless I've > done something wrong > by selecting the regions I want to plot with departements[-exclude]. I think here the problemis not using exact=T in the call to map(), see below. > I also get an error > when I try to use map.text to label a map with only the regions I'm > selecting. > > > departements <- map('france',namesonly=T, plot=FALSE) > > # returns a vector of names of regions > > > > exclude <- c(47, #Alpes-Maritimes > + 66, # Haute-Savoie > + 76, # Savoie > + 95, # Territore-de-Belfort > + 109, 110, 111, # Var: Iles d'Hyeres > + 49, 53, 54, 55, # Moribhan: Isles > + 62, 64,# Vendee: Isles > + 72, 75 # Charente-Maritime: Isles > + ) > > > > depts <- departements[-exclude] > > gfrance <-map('france', regions=depts) > > labels <- (as.character(1:length(departements)))[-exclude] > > gfrance <-map.text('france', regions=depts, add=FALSE, labels=labels) > Error in map.text("france", regions = depts, add = FALSE, labels = labels) : > map object must have polygons (fill=TRUE) > That error message is issued when regions= specifies less than the whole database, in which case the alternate "list of 'x', 'y', and 'names' obtained from a previous call to 'map'" is required as a first parameter. So to do what you want (notwithstanding the next problem you raise), try the following as a demonstration: gfrance <- map('france', regions=depts, exact=T, fill=T, plot=F) map.text(gfrance, regions=depts, labels=labels) map('france', regions=departements[exclude], fill=T, col=1, add=T) > Another problem, potentially more difficult for mapping data on the map > of France is that > the "departements" are actually just the polygons in the map, > arbitrarily numbered from > east to west, and from north to south --- they don't correspond to the > 'official' administrative > region numbers. As well, the departement names don't always match > exactly (ignoring > accents, e.g., Val-d'Oise vs. Val-Doise) so it would be another > challenge to plot my > historical data on the map of France. > Well, maps is a source package! [:-)]. You are most welcome to modify the source files maps/src/france.{gon,line,name} to reorder the polygons (and correct errors in the names). If the relationship between those 3 files is not obvious, contact me for further details. Also, I am happy to fold your changes back into the original maps package. Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] adjusting the map of France to 1830
At 16:29 18/11/2004, Michael Friendly wrote: > I'm doing some analyses of historical data from France in 1830 on 'moral > statistics' that I'd like to > show on a map. I've done most of my analyses in SAS, but a few things > would work better in R. > To do this, I have to adjust the modern map, > > library(maps) > map('france') > > to adjust for changes in departments (86 in 1830, to 97 now). I've read > the documentation > for the maps and maptools package, but there seems to be no functions to > allow this, and > I can't find information on the exact structure of map datasets, but I > understand them to > be delimited lists of polygon coordinates. > > In SAS, all maps have (one or more) ID variables representing the > geographical region, > and there is also a proc gremove that can remove internal boundaries > inside the polygons > for regions with the same ID. Is there some way I can do this in R? > Unfortunately not with the current implementation of several of the 'extra' databases in the mapdata package. The map() function does have the interior=FALSE option, which would normally do what you want, but only when the data has been manipulated to allow it. Currently this option is only useful with the world and usa maps (and their derivatives, such as world2 and state). Currently every department is a complete polygon, and so every interior line segment occurs twice in the data. What has to happen to the data is for it to be split up into non-overlapping line segments, and each polygon reconstructed from a list of these line segments (with direction being important). If you are prepared to perform this somewhat tedious process, I am happy to assist you with further details. However even with the interior= option functioning, it would still not be easy to produce the map you would require, since you would have to build it up from many components (namely each of the 'combined' departments, plus 'all the rest'). HTH, Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] Seg. faults in mapthin (in package maps!)
> From: "Kevin Bartz" <[EMAIL PROTECTED]> > Date: Thu, 1 Jul 2004 10:18:36 -0700 > > Dr. Ripley, you are the master. That fix worked like a charm! All the way to > 50, with no problems. Thanks again, > Hear, hear! > From: Prof Brian Ripley [mailto:[EMAIL PROTECTED] > Sent: Thursday, July 01, 2004 10:07 AM > > There are a few other quirks, and I would check some of the other uses of > `long' in the C sources. > OK, I've changed most of the 'long's to 'int's in the next version (maps_2.0-22), soon to appear on CRAN. It passes R CMD check on my 32-bit systems, but please test it thoroughly on a 64-bit system. Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] Seg. faults in mapthin (in package maps!)
> From: "Kevin Bartz" <[EMAIL PROTECTED]> > Date: Thu, 1 Jul 2004 10:18:36 -0700 > > Dr. Ripley, you are the master. That fix worked like a charm! All the way to > 50, with no problems. Thanks again, > And all while I was asleep! I'll submit a corrected maps package in a couple of weeks, when I return from China (no guaranteed connection where I'll be). Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Greater than 1 or less than 1?
> Date: Sun, 20 Jun 2004 18:32:46 -0400 > From: "Zhen Chen" <[EMAIL PROTECTED]> > > I have problem evaluating the expression h = exp(x)/(exp(exp(x))-1) for > large negative x. This expression is actually the probability > that y = 1 when y is a Poisson random variable truncated at 0, hence > must satisfy 0 <= h <= 1. However, when > x < -18, I may get an h value that is larger than 1 while the true > value should be a quantity that is smaller than but very close to 1. > For example when x = -19, h = 1.031174. I tried to use > different form of h and none of > them give me an h value that is less than or equal to 1. I also tried > to find patterns, but discovered none: for some x < -18, h is below 1; > for others > h is greater than 1. Is there any trick that enables me to obtain > theoretically correct results from this expression? > exp(x)/(exp(exp(x)) - 1) is less than 1 for me for x = -18 with R-1.9.0 on a Solaris 9 system, a NetBSD 2.0C system and a Windows XP system. You could try exp(x)/expm1(exp(x)), but if the longer expression doesn't work, then this may not either, though it does give slightly different results on my systems. E.g.: > x = -19 > exp(x)/(exp(exp(x)) - 1) - 1 [1] -2.047911e-09 > exp(x)/(expm1(exp(x))) - 1 [1] -2.801398e-09 log(1 + x) and exp(x) - 1 are often treated specially in math libraries for when x << 0. Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] R on Tru64 OSF 5.1
> From: =?iso-8859-1?q?J=20Swinton?= <[EMAIL PROTECTED]> > > Has anyone achieved a current successful build of R > 1.8 or 1.9 or earlier for HPUX Tru 64 OSF 5.1? > > There is no binary version of R for HPUX Tru 64 OSF > 5.1. The R admin manual mentions that the native make > fails on "Alpha/OSF (aka Tru 64)" and gnu make must be > used instead. There are problems reported with > building R versions around 1.4 from 2002 on R-devel, > some of which were fixed with using gmake instead of > make, some of which were fixed by changes to R > sources, but it is unclear if current R compiles on > the current OSF 5.1, which is what I need to know > before trying to get our local IS admins to try again > after their first three failures. > > Jonathan > It depends what you mean by "successful build". I succeeded in getting R-1.8.1 to "./configure" and "make" on a system: Compaq Tru64 UNIX V5.1A (Rev. 1885); Wed Nov 20 00:37:24 NZDT 2002 However it does not pass "make check", with the following message: 496306:/vol/R/src/osf1/R-1.8.1/bin/R.bin: /sbin/loader: Fatal Error: call to unresolved symbol from /vol/R/src/osf1/R-1.8.1/modules/lapack.so (pc=0x3ffbfde3308) during reg-tests-1. This is possibly a consequence of the following (Info) warnings issued at make time: f77 -fpe3 -g -c dlapack0.f -o dlapack0.lo f90: Info: dlapack0.f, line 1721: This directive is not supported in this platform. CDEC$ NOVECTOR ---^ f77 -fpe3 -g -c dlapack1.f -o dlapack1.lo f77 -fpe3 -g -c dlapack2.f -o dlapack2.lo f77 -fpe3 -g -c dlapack3.f -o dlapack3.lo f90: Info: dlapack3.f, line 21481: This directive is not supported in this platform. CDEC$ NOVECTOR ---^ f90: Info: dlapack3.f, line 21495: This directive is not supported in this platform. CDEC$ NOVECTOR --^ f90: Info: dlapack3.f, line 21522: This directive is not supported in this platform. CDEC$ NOVECTOR ---^ The weird thing is that the problem occurs during the following code fragment: ## predict.arima0 needed a matrix newxreg: Roger Koenker, 2001-09-27 u <- rnorm(120) s <- 1:120 y <- 0.3*s + 5*filter(u, c(.95,-.1), "recursive", init=rnorm(2)) fit0 <- arima0(y,order=c(2,0,0), xreg=s) which runs perfectly OK if pasted into an interactive session. I did spend some time trying to isolate/reproduce the problem without success. Unfortunately, this is the only Alpha/OSF1 system we have and it is not used much (particularly not for R), so the installation has not been well tested for 'real use'. Let me know if you want the necessary configuration changes. Hope this helps, Ray Brownrigg <[EMAIL PROTECTED]> http://www.mcs.vuw.ac.nz/~ray __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] writing text on graphics' window
> From: Duncan Murdoch <[EMAIL PROTECTED]> > Date: Sun, 21 Mar 2004 06:58:03 -0500 > > On Sun, 21 Mar 2004 11:07:01 +, you [EMAIL PROTECTED] wrote: > > >Does anyone know of a method for writing text to the graphics window, > >where there is *no* plot? Basically, I have developed a 'significance > >test' and I would like the output on the graphics window to say > >something about the input parameters and the stats of the significance test. > > You need to make sure a graphics device is active and establish a > coordinate system there. The easiest way to do that is to make a call > to plot() with everything turned off: > > plot(0:100,0:100,type='n',axes=FALSE,xlab="",ylab="") > > You may also want to reduce the "margins" if you want your output to > take up the full frame, e.g. > > oldmargins <- par(mar=c(0,0,0,0)) > plot(0:100,0:100,type='n',axes=FALSE,xlab="",ylab="") > An easier way to activate a graphics device and establish a coordinate system is to call plot.new(). Try: > plot.new() > text(0, 0, "ABC") > par("usr") [1] -0.04 1.04 -0.04 1.04 Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] R equiv to proc gremove in maps package
> Is there an R equivalent to SAS's proc gremove? You would use this procedure to > combine the units on an existing map, for example to build a map of Metropolitan > Statistical Areas (MSAs) from the [US] counties dataset where the internal boundries > surround the MSAs (which are groups of counties) rather than the individual > counties. I can imagine the mechanism would be to find and erase the boundries > within each group. Is there a way to do that in R? I looked at the S reference on > building geographical databases, and I could probably fumble through the C code to > build the binary files of lines and polygons, but have no idea how to implement the > shell program to extract the datapoints. Thanks. > It's already built in to the maps package. Try e.g.: map("county", c("arizona,gila", "arizona,maricopa", "arizona,yavapai"), interior=F) [I wouldn't have believed that it worked until I tried!] Hope this helps, Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Map of British Colonial America 1775
> I've reviewed the "map" documentation I've found so far with > S-Plus 6.2 and R including "www.r-project.org" -> search -> "R site > search" -> "world map" and Becker & Wilkes (1993) "Maps in S" cited in > the "map" documentation. I've noticed that, for example, > 'library(maps); map("world", "UK:Bermuda")' produces an 83-point > polygon sketch in S-Plus 6.2 but only a degenerate 3-point polygon (with > the third point = the first) in R 1.8.1. > You need the mapdata package. The default "world" database in S-Plus is the high-resolution one (with also "world.thin"). In R the high resolution one is "worldHires" from mapdata (mainly so that the standard maps package is a manageable size). Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] drawing filled countries according to data using map('world')?
> My question then is, how to compute a similar procedure using the 'world' > database. Specifically, how can I access the country names in the 'world' > database to accomplish the translation to the country names in my dataset? > Is there any way to unpack the 'world' database to do the matching in an > external program? And does anybody now of other (more recent) world maps > that I could use? > The short answer is that the file world.N in .../library/maps/mapdata/ contains the mapping from polygon numbers to names. I think the only other thing you need to know is that the colours are allocated sequentially in increasing numerical order of polygon number. Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] ISO country codes available in R's {map}?
> Does any one know wether it is possible to get the Country codes from the > Map package - I've not found anything in the docs and only seen a post re. > how to get country names. I'd like to programmatically marry some data > indexed by ISO country code to the map package contry names. If not then > I'll need to hunt for a listing and try and marry that to the {map} > names. > No part of the maps package knows anything about ISO country codes. I guess if there was a database of ISO country codes somewhere, it wouldn't be too difficult to provide some add-on functionality to merge the internal polygon numbers with country code numbers. The hardest part would be what to do with mismatches. Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] Distance and Aggregate Data - Again...
[EMAIL PROTECTED] wrote: > I appreciate the help I've been given so far. The issue I face is > that the data I'm working with has 53000 rows, so in calculating > distance, finding all recids that fall within 2km and summing the > population, etc. - a) takes too long and b) have no sense of progress. > Well, there are ways to speed this up. > Below is a loop that reads each recid one at a time, calculates the > distance and identifies the recids that fall within 2 km. It iterates > through all records successfully. > But you don't need to subset d[, 2:3] every time, e.g. Also, my experience is that writing a single line appended to a file every time around the loop is very inefficient. > Where I'm stuck is how to get the sum of population and dwellings and > the mean age for the records that are selected. Also, the desired > output should have the following fields: recid, sum(pop), sum(dwell), > mean(age). I don't know how to write only those fields out to the > file. > Well, that part is easy, essentially, you have to save the indices, then extract the relevant rows of pop, dwell and age. So I would modify your code as follows: library(fields) d <- as.matrix( read.csv("filein.csv") ) lonlat2 <- d[,2:3] results <- matrix(0, nrow=nrow(d), ncol=4) for(i in 1:nrow(d)){ lonlat1 <- d[i, 2:3] whichdist <- which(rdist.earth(t(as.matrix(lonlat1)), as.matrix(lonlat2), miles=F) < 2) distval <- d[, 1][whichdist] sumpop <- sum(data[, "pop"][whichdist]) sumdwell <- sum(data[, "dwell"][whichdist]) meanage <- mean(data[, "age"][whichdist]) results[i, ] <- c(d[i, "recid"], sumpop, sumdwell, meanage) } write(t(results), file="C:\\outfile.out", ncol=ncol(results)) Then there is a trick you can use to speed up the process. Essentially you reduce the 'inner loop' which is inside rdist.earth(). This is achieved by initially computing a single distance vector with reference to a fixed point [(0, 0) seems reasonable since it is in the Atlantic Ocean]. Then you select a subset of your points that are the same distance is your circle centre from the fixed point plus or minus the radius you want (and a bit of tolerance). This will generate an annulus of points rather than a circle, but in particular all the points within the circle of interest will also be within this annulus. Then you apply rdist.earth() to find the distance of each of these points from your circle centre. This is what I have used to process a 52000 length matrix in about 40 minutes (and ~50MB memory): DistAgg <- function(data, dist=2) { results <- matrix(0, nrow=nrow(data), ncol=4) colnames(results) <- c("recid", "sumpop", "sumdwell", "meanage") lonlat2 <- data[, 2:3] basedist <- rdist.earth(t(as.matrix(c(0, 0))), as.matrix(lonlat2)) for(i in 1:nrow(data)){ lonlat1 <- data[i, 2:3] approxval <- which(abs(basedist - basedist[i]) < dist*1.001) if (length(approxval) > 1) { whichval <- approxval[which(rdist.earth(t(as.matrix(lonlat1)), as.matrix(lonlat2[approxval, ]), miles=F) < dist)] } else { whichval <- approxval } sumpop <- sum(data[, "pop"][whichval]) sumdwell <- sum(data[, "dwell"][whichval]) meanage <- mean(data[, "age"][whichval]) results[i, ] <- c(data[i, "recid"], sumpop, sumdwell, meanage) } write(t(results), file="C:\\outfile.out", ncol=ncol(results)) } The data I used was based on yours, but randomised latitude and longitude (in restricted ranges). Hope this helps, Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] parse error in GLMM function
"CENDOYA, Gabriela" <[EMAIL PROTECTED]> writes: > I'm trying to use the function GLMM from lme4 package, (R-1.8.1, Windows > 98),and I get the following error: > > pd5 = GLMM(nplant~sitio+ > + fert+ > + remo+ > + sitio:fert+ > + remo:sitio+ > + remo:fert+ > + remo:fert:sitio > + data=datos, > + family=binomial, > + random=~repe:sitio) > Error in parse(file, n, text, prompt) : parse error Isn't the problem that you really do have a syntax error? Note there is no comma before the data=. Ray __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] large fonts on plots
On 11 Feb 2004 at 13:30, Piet van Remortel wrote: > Hi all, > > I need to enlarge te fonts used oo R-plots (plots, histograms, ...) in > labels and titles etc. Nobody seems to have mentioned another way of achieving this, which is to scale the original plot to suit the intended final size in the document. E.g. if you use: > postscript(file="plot1.ps", width=4, height=4) then the 4"x4" figure will 'automatically' have relatively larger text and plot symbols. This works well with x11() which is subsequently printed, and pdf() also. HTH, Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Re: [R] filled maps
> Date: Tue, 3 Feb 2004 11:34:06 +0100 (CET) > From: Roger Bivand <[EMAIL PROTECTED]> > > On Tue, 3 Feb 2004, Janus Larsen wrote: > > > Hi R-Help, > > > > I would like to make filled contour maps of ocean data overlaid by > > costlines from the map package. > > I can draw the filled contours and the coastlines om the same plot, but > > the filled contour also covers part of the land. To get rid of that I > > tried to draw a filled coastline map on top of the filled contour, but > > the filled map only draws the closed contours - so most of the land is > > missing. > > Example: > > map("worldHires",xlim=c(0,15),ylim=c(50,60)) #Draw relevant region > > (North Sea and Denmark waters) > > map("worldHires",xlim=c(0,15),ylim=c(50,60),fill=TRUE) # This only draws > > Denmark and Holland (Sweden, uk, Germany etc. disappears because they > > are not closed polygons). > > Looks as though you can use the regions= argument: > > > res <- map("worldHires",xlim=c(0,15),ylim=c(50,60), plot=FALSE, > + namesonly=TRUE) > > map("worldHires",xlim=c(0,15),ylim=c(50,60)) > > map("worldHires",regions=res, fill=TRUE, add=TRUE) > Thanks Roger, I hadn't got around to figuring out that workaround. Certainly the mapgetg code explicitly excludes incomplete polygons if fill=TRUE, but since the polygon() function is happy to clip what it is filling, your workaround takes advantage of that. I'll have a look at getting mapgetg to do this itself - it would have to always provide complete polygons if fill=TRUE instead of just 'visible' polylines. > fixes it for me. Curiously, the longitudes add 2.5 degrees on each side. > Actually this is to preserve aspect ratio. An earlier version of the maps package adjusted the shape of the figure to conform to both aspect ratio and imposed limits. I didn't deliberately change this, but I think I know why it has happened, so I'll look at this also. Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
RE: [R] beginner programming question
> From: "Gabor Grothendieck" <[EMAIL PROTECTED]> > Date: Wed, 17 Dec 2003 15:02:49 -0500 (EST) > > Define function f to take a vector as input representing > a single input row. f should (1) transform this to a vector > representing the required row of output or else (2) produce > NULL if no row is to be output for that input row. > > Then use this code where z is your input matrix: > > t( matrix( unlist( apply( z, 1, f ) ), 2) ) > But as has been pointed out recently, apply really is still just a for loop. > > From: Adrian Dusa <[EMAIL PROTECTED]> > > Date: Wed, 17 Dec 2003 21:28:05 +0200 > > > > I have a (rather theoretical) programming problem for which I have found > > a solution, but I feel it is a rather poor one. I wonder if there's some > > other (more clever) solution, using (maybe?) vectorization or > > subscripting. Here is a subscripting solution, where (for consistency with above) z is your data [from read.table(filename, header=T)]: > z rel1 rel2 rel3 age0 age1 age2 age3 sex0 sex1 sex2 sex3 113 NA 25 232 NA121 NA 2413 35 67 34 102212 3144 39 40 59 601221 44 NA NA 45 70 NA NA22 NA NA > res <- matrix(NA, nrow=length(z[, 1]), ncol=2, dimnames=list(rownames=rownames(z), colnames=c("ageh", "agew"))) > w <- w0 <- w1 <- w2 <- which(z[, c("rel1", "rel2", "rel3")] == 1, T) # find spouse entries > w0[, 2] <- z[, "sex0"][w[, 1]]# indices for respondent's age > w1[, 2] <- 3 - w0[, 2]# indices for spouse's age > w2[, 2] <- 4 + w[, 2] # indices of spouse's age > res[w0] <- z[, "age0"][w[, 1]]# set respondent's age > res[w1] <- z[w2] # set spouse's age > res colnames rownames ageh agew 1 25 23 2 34 35 3 39 40 4 NA NA > Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
Re: [R] Recoding problem
> I have the following variables, all of which are logicals > > fmar15 fcoc15 fher15fcrk15fidu15 > > what I would like is a variable drug15 which equals > idu if fidu15 is T; crk if fidu15 is F but fcrk is T, her if fher15 is > T but fcrk15 and fidu15 are F and so on > > What's the best way to do this? > I don't know about the best way, but if I understand your question, the following achieves what you want: > tab <- cbind(fidu15, fcrk15, fher15, fcoc15, fmar15) > substring(colnames(tab), 2, 4)[apply(tab, 1, match, x = T)] HTH Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
Re: [R] ?for
> I thought I was getting to grips with it, until last reply from Prof Brian > Ripley's. At least in my Windows and Linux boxes Actually it was my reply that suggested ?paste("help") > ? paste("help") > opens the help for help (!?) Not on my systems! [Solaris, NetBSD, Windows 98]. I believe it should open the help for paste(). The point is that it requires a syntactically correct statement after the ? so that the command parses correctly (with ? as a unary operator), but the result is to provide the help for the first syntactic element of the parsed expression, *without* evaluating the expression. Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
Re: [R] file not found?
Thomas W Blackwell <[EMAIL PROTECTED]> wrote: > I recall a message on the [r-pkgs] list from Ray Brownrigg > on November 1 this year that may be relevant. At that time > he was announcing the availability of maps_2.0-8.zip and > some other packages. > > Therefore, 2.0-10 must be *very* recent and might still be > going through the CRAN propagation process. (I don't know > anything about whther it is or not.) I would wait a couple > of days and try again, if I were you. > Well, actually, 2.0-10 in source form is on CRAN since last week, but for some reason only in the 1.7/ directory in the Windows area. I just this morning noticed this and sent a message to the Windows CRAN maintainer. In the meantime, 2.0-11 has been uploaded, probably to appear in a few days (no noticeable changes). Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
Re: [R] Problem with parser and if/else
Prof Brian Ripley <[EMAIL PROTECTED]> > On Fri, 14 Nov 2003, Paul Lemmens wrote: > > > I'm trying to grasp this: if you're saying (or are you saying) that the > > only way to have if() know that an else will be present is to put it on the > > same line as the closing curly brace of the if() (compound) statement. But > > if I look at some code from, e.g., aov and lm, I see plenty violations of > > that rule. > > The actual rule is given in my reference (the one that Ben Bolker did not > bother to look up) earlier in this thread. > > You need to ensure that the code is syntactically incomplete when `else' > is encountered. That will always be true inside a braced expression such > as the bodies of the functions you quote. But at top-level, you do need > to write > > if(condition) something else something_else > > or > > if(condition) { > } else { > } > > since > > if(condition) { > } > else { > } > > fails the test. > BUT, just to make it even clearer, what does succeed interactively is (what is effectively the same as a function body): { if(condition) { "TRUE" } else { "FALSE" } } because the parser has to find the final closing brace before the syntactic item is complete. Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
Re: [R] stop further sourcing of an R file
> Thankyou, but i *donot* want the R session to quit. > Someone suggested writing functions with return kept at > the locations. > I wanted to ask if there was another way. > Just have a single } as your last line, and put: if (FALSE) { just after where you want to stop. Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
Re: [R] Test for new page
Marc Schwartz <[EMAIL PROTECTED]> wrote: > On Tue, 2003-11-11 at 16:49, Paul Murrell wrote: > > Hi > > > > > > Aleksey Naumov wrote: > > > Dear R experts, > > > > > > I am writing a multi-page PDF file and would like to put a header on each > > > page. Is there a way to test a graphic device to see if a new page is started > > > (so that I know when to write the header)? > > > > > > Sorry. Not that I can think of. > > > > Paul > > > > > > > I could simply count the plots made (each page has the same number of plots), > > > but wanted to see if a more general solution is available. > > > I was trying to think of a way but could not either. > It's not clear what you want as the "header", but if it is the same on each page, then using mtext("Header", outer=TRUE) (after a suitable par(oma=c(0, 0, 1, 0))), you can just write the header after every plot. Ray __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
Re: [R] map does not display maps, MacOSX
Philippe Glaziou <[EMAIL PROTECTED]> wrote: > > I installed the maps and mapdata libraries on my R-1.8.0 on > MacOSX 10.2.8 (jaguar on a powerbook G4), and failed to make the > map function work properly: > > > R : Copyright 2003, The R Development Core Team > Version 1.8.0 (2003-10-08) > [...] > > > library(maps) > > map() > Error in file(file, "r") : unable to open connection > In addition: Warning message: > cannot open file `/Users/glaziou/Library/R/maps/mapdata//world.N' > > map('usa') > Error in file(file, "r") : unable to open connection > In addition: Warning message: > cannot open file `/Users/glaziou/Library/R/maps/mapdata//usa.N' > > system('ls -l /Users/glaziou/Library/R/maps/mapdata') > total 1796 > -rw-r--r-- 1 root staff 143902 Oct 14 11:30 county.G > -rw-r--r-- 1 root staff 690260 Oct 14 11:30 county.L > -rw-r--r-- 1 root staff 618 Oct 14 11:30 nz.G > -rw-r--r-- 1 root staff 13040 Oct 14 11:30 nz.L > -rw-r--r-- 1 root staff2642 Oct 14 11:30 state.G > -rw-r--r-- 1 root staff 96892 Oct 14 11:30 state.L > -rw-r--r-- 1 root staff 282 Oct 14 11:30 usa.G > -rw-r--r-- 1 root staff 58232 Oct 14 11:30 usa.L > -rw-r--r-- 1 root staff 74434 Oct 14 11:30 world.G > -rw-r--r-- 1 root staff 295152 Oct 14 11:30 world.L > -rw-r--r-- 1 root staff 74434 Oct 14 11:30 world2.G > -rw-r--r-- 1 root staff 295152 Oct 14 11:30 world2.L > -rw-r--r-- 1 root staff 54832 Oct 14 11:30 world2.N > Note the actual file requested (usa.N) doesn't actually exist! > : > > Any hint appreciated, > I had this sort of problem on a Windows system, and had to add an extra command into Makefile.win. [Is there recognised such a thing as Makefile.mac? I seem to recall there used to be, but it is not mentioned in the latest R-exts.] What you need to do is add the line: $(CP) ${*}.n ../inst/mapdata/${*}.N # need this here for Mac as the second line of the .gon.g: target in maps/src/Makefile. [note the first character is a tab] I'll fix this (and the //) in the next version. Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
Re: [R] USA map
> Ivar Herfindal <[EMAIL PROTECTED]> wrote: > > On Wed, 5 Nov 2003 11:28:42 +0100 (MET), Mathieu Ros > <[EMAIL PROTECTED]> wrote: > > >>>>>> "k" == kjetil <[EMAIL PROTECTED]> disait: > > > > > > k> I also tried > > > > k>> map("worldHires","sweden") > > k>> map("worldHires","denmark") # which comes out very small since it > > k> # includes the Faroe k> > > # islands properly faraway > > > > and, just to know, how would you do to plot *only* continental > > denmark? The same applies for france, UK, ... > > > Hello > > One simple way of doing it is to specify the xlim and ylim of your map, > e.g. library(maps) > map('world', 'Norway', xlim=c(5, 33), ylim=c(55, 75)) > But the 'best' way is RTFM! map("world", "Norway", exact=T) Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
Re: [R] USA map
> From [EMAIL PROTECTED] Wed Nov 5 11:31:52 2003 > > but I cannot find any function to find the avalable region names in > a map database > They are in the .N file for the appropriate map. Try: library(maps) data(worldMapEnv) file.show(paste(Sys.getenv(worldMapEnv), "world.N", sep="/")) These ASCII files are in the form: region[:subregion] polygonNumber Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
Re: [R] USA map
Prof Brian Ripley <[EMAIL PROTECTED]> wrote: > > On Tue, 4 Nov 2003 [EMAIL PROTECTED] wrote: > > > On 3 Nov 2003 at 15:04, [EMAIL PROTECTED] wrote: > > > > > > > > I believe in R you need > > > > > > library(maps) > > > > > > before using usa(). Of course, the library has to be installed first. > > > > > > > Yes, I new that. When I did help.search("usa") I HAD maps > > attached. It dosen't have any function usa(). R is more international > > minded: > > > > > library(maps) > > > usa() > > Error: couldn't find function "usa" > > > map('usa') > > > > > map('china') > > map('bolivia') > > > > and so on also works! > > You need library(mapdata) for those. > Except the Bolivia map must be a local one (from [EMAIL PROTECTED]); it is not in library(mapdata). Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
Re: [R] USA map
> > I don't know the usa() function in S but usa() in R gives a map of the US > > Where did you find that function? In rw1080 on windows XP I get: > Some earlier versions of the maps package (which was then available for "Unix" only), did have individual functions for each country in the world database. In the current version, you need to type: library(maps) map("usa") #or: map("world2", "USA") Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
[R] [R-pkgs] R maps, mapdata and mapproj for Windows
Nearly a month ago I announced the availability of a Windows port of the maps for S package but for various reasons the Windows binaries did not appear on CRAN until yesterday. Now I am pleased to announce that the three packages maps_2.0-8 (the base maps package with low-resolution databases), mapdata_2.0-6 (high-resolution and other databases) and mapproj_1.1 (projection add-on for maps, maintained by Tom Minka) are all now available on CRAN both in source and Windows binary forms. I apologise for the delays. Ray Brownrigg ___ R-packages mailing list [EMAIL PROTECTED] https://www.stat.math.ethz.ch/mailman/listinfo/r-packages __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
[R] [R-pkgs] R maps now for Windows
The original S maps package, previously only available in R on Unix systems is now available in R for Windows. maps_2.0-2 and mapdata_2.0-2 (see below) in both source and Windows binary form (for rw1071) have been uploaded to CRAN, hopefully soon to appear. The key aspect of this package which differentiates it from other GIS related R packages under development is that it is standalone within R, i.e. it does not rely on the presence of any other software. Features are: 1) Built-in databases for usa, state (states of the USA), county (counties of the USA), world, world2 (Pacific-centric world) and nz. 2) Databases supplied in source form with accompanying software to generate the binary forms, allowing users to edit the data or add new databases (basic format is latitude/longitude pairs). 3) A database of world cities (generally >~100,000 population). 4) A companion package called mapdata provides high-resolution versions of world, world2 and nz, as well as china (provinces) and rivers. 5) A companion package called mapproj provides a variety of map projections (from Doug McIlroy). 6) While being backwards-compatible with the original 'S maps library', it also has many new features, including: - displaying region names on a map, - finding the map region containing a given point (point-in-polygon), - computing areas of map regions, and - interpolating data aggregated over map regions (for statistical maps). The maps and mapdata packages are maintained by Ray Brownrigg with generous assistance from Tom Minka. The mapproj package is available at http://www.stat.cmu.edu/~minka/software/maps/ Ray Brownrigg Tom Minka October 2003 ___ R-packages mailing list [EMAIL PROTECTED] https://www.stat.math.ethz.ch/mailman/listinfo/r-packages __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
Re: [R] trouble with maps
[EMAIL PROTECTED] wrote: > On Wed, 23 Jul 2003, Deborah Swayne wrote: > > > Has anyone else seen this behavior from the "maps" package? > > > > map('state', fill=TRUE) > > > > results in a lively mix of overlapping polygons inside a > > map of the US, but they have no obvious relationship to > > state boundaries. (See attached jpeg.) > Ah, this was a 'known problem' :-) first noticed in January this year, by Ott Toomet <[EMAIL PROTECTED]>. I suspect it is a feature of the port to R which never worked properly. The solution is to reverse the sign of every line number making up the polygons in the .gon file (or alternately, reverse the order of the line numbers there). A 'fixed' maps package (Unix only) is available at: ftp://ftp.mcs.vuw.ac.nz/pub/statistics/map/maps_1.1-2.tar.gz Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
Re: [R] maps library for R?
> Does it already exist a library in R to draw maps (something like a Generic > Mapping Toolbox, http://gmt.soest.hawaii.edu/ port). > > I've seen in an old R-help that Ross Ihaka once tried to port the S-plus map > library (http://maths.newcastle.edu.au/~rking/R/help/99b/0832.html ). > Anybody know if this package is available somewhere or if there is somebody > developing a Mapping package for R. I'll might then try to help. There are two, but AFAIK neither of them will work on Windows. My port of the basic S maps package (no projections) has recently been modified to provide a small (700KB) 'base' package with limited maps data (just usa, state, counties, nz and world.thin), plus an add-on package called mapdata (24MB) with the high-resolution maps. [This split has not yet been fully tested, and requires maps to be *installed* before mapdata, and mapdata must be installed in the same library directory as maps.] These are available from ftp://ftp.mcs.vuw.ac.nz/pub/statistics/map/map*_1.1-0* The other solution in development is called Rmap from Barry Rowlingson (http://www.maths.lancs.ac.uk/Software/Rmap/), but as I understand it, this requires quite a lot of "third-party" stuff to be installed first. Ray Brownrigg __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
Re: [R] Combining the components of a character vector
>Suppose I have a character vector. > > x <- c("Bob", "loves", "Sally") > > I want to combine it into a single string: "Bob loves Sally" . > paste(x) yields: > paste(x) > [1] "Bob" "loves" "Sally" > > Is there a more natural (no loop) way to do this in R? > RTFM: > paste(x, collapse=" ") [1] "Bob loves Sally" Ray __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
Re: [R] inserting elements in a list
> From: Peter Wolf <[EMAIL PROTECTED]> > > "insert.values" <- function(x,pos.insert,x.insert){ > : Which when reduced to its bare minimum to solve the stated problem, becomes: threes <- which(a==3) a <- c(a, rep(7, length(threes)))[order(c(seq(a), threes))] This is by far the overall fastest solution so far proposed, and certainly is reasonably easy to understand. However an even faster solution is based on rep(): N <- which(threes <- a == 3) # which ones a <- rep(a, times = 1 + threes)# duplicate them a[N + 1:length(N)] <- 7# replace the duplicates Relative timing for vector lengths n=1000 and n=1 with 1000 repeats (on a 1.7GHz P4 running R-1.6.1 under NetBSD) is as follows: a <- sample(1:9, n, TRUE) n=1000 n=1 lapply (Stephen Upton and Peter Dalgaard) 9.2s 117.8s recursion (John Fox) [options(expressions=n)]59.9s stack overflow offset <- (Thomas Lumley) 1.4s 12.7s offset <- (Peter Dalgaard)1.3s 10.9s order() (Peter Wolf) 1.1s8.7s rep() (me) 0.9s7.6s Hope this helps, Ray Brownrigg <[EMAIL PROTECTED]> http://www.mcs.vuw.ac.nz/~ray __ [EMAIL PROTECTED] mailing list http://www.stat.math.ethz.ch/mailman/listinfo/r-help
Re: [R] Interpolation
> I have two vectors and connect the points by line segments: > > x <- c(1990,1994,1995,1997) > y <- c(30,40,80,20) > plot(x,y,type="l") > > I want to get the values of y's on these lines at missing x's, > i.e., at 1991,1992,1993,and 1996. Is there a simple way to do this? > Easy: approx(x, y, c(1991:1993, 1996))$y while approx(x, y, 1990:1997)$y will give you all the y values. Hope this helps, Ray Brownrigg <[EMAIL PROTECTED]> http://www.mcs.vuw.ac.nz/~ray __ [EMAIL PROTECTED] mailing list http://www.stat.math.ethz.ch/mailman/listinfo/r-help