gmane.comp.lang.r.geo

Re: raster package options maxmemory and chunksize

Robert J. Hijmans — 2012-05-18T01:50:58

Thanks Tim,

That is good to know. The settings are somewhat conservative to avoid
virtual memory use on macs (which slows things down a lot on that
platform). Still, even the default settings should beat scanning ascii
files, so there might be room for improvement. One reason might be that the
raster::predict function tries to accommodate many different
models/situations, so there is more overhead.

I do not know if the native file format helps here. My guess is that it
does not, but that it might if you used a single multi-layer file
(RasterBrick) in the native format.

In your case, the most limiting factor is not the file I/O but the
randomForest prediction. With 8 cores I see a reduction in run time of up
to 50%. If you are doing regression (classification may not work) all you
need to do is run

beginCluster()

Robert



On Thu, May 17, 2012 at 5:05 AM, Tim Howard <tghoward< at >gw.dec.state.ny.us>wrote:


[[alternative HTML version deleted]]

Plotting cross-correlograms

Michael Hewson — 2012-05-17T23:04:09

G'day List,

I'm using package "ncf" to compute and plot cross-correlograms. The documentation for plot.correl notes the syntax is plot.correlg(x, ...) with "x" being the item of class correlog and "..."being "other arguments".

I'm interested in finding out what I can put into "other arguments" please. Are the permissible arguments from another R package?

I'm actually trying to plot two cross-correlograms on the same graphic and combined axes. I'm guessing I might have to construct my own class correlog plot routine from other packages to do that?

With thanks,

Michael Hewson PhD candidate
Climate Research Group<http://www.gpem.uq.edu.au/crg> | Centre for Spatial Environmental Research<http://www.gpem.uq.edu.au/cser>
The University of Queensland | Brisbane Q 4072 | Australia





[[alternative HTML version deleted]]

how to determine varying spatial correlation scale

m p — 2012-05-17T19:45:33

Hello,
I used "correlogram" from "spatial" package to determine correlation scale
for my 2d data (x,y,value)  but just looking with bare eye it seems that
the correlation scale varies over the domain.
How can I best determine regional correlation scales? Can someone suggest
what would the best way to handle that problem?
Thanks,
Mark

[[alternative HTML version deleted]]

Re: Raster getData() problem

Robert J. Hijmans — 2012-05-17T19:08:34

Raúl,

You are not doing anything wrong. The function does not work for the few
countries for which the rasters have been split in multiple files. It is
has been fixed for the next release or 'raster'.

Thanks, Robert

On Tue, May 15, 2012 at 12:58 AM, Raúl Jiménez Rosenberg <
raul.jimenez.rosenberg< at >gmail.com> wrote:


[[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
R-sig-Geo< at >r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: Calculating area of region with a raster

Robert J. Hijmans — 2012-05-17T19:04:17

Rhys,

This returns the approximate area of each cell, which varies with latitude
as you use longitude/latitude data.

a <- area(ra, na.rm=TRUE)

to sum the area by zone (RasterLayer 'world' in your case), you can do

zonalarea <- zonal(a, world)



For a Raster object with a _planar_ (not Longitude/Latitude) coordinate
reference system, the area function is not relevant as each cell has the
same size, and you could do

f <- freq(world)

and, to go from cellcounts to areas:

f[,2] <- f[,2] * prod(res(world))

This is what Tony suggested for a RasterLayer with presences/absence data.
The memory-safe version of his approach is:

icearea <- cellStats( raster.Ice , 'sum', na.rm=TRUE)  * prod(res(
raster.Ice ))



Robert


On Wed, May 16, 2012 at 3:56 PM, Anthony Fischbach <afischbach< at >usgs.gov>wrote:



Hi all,

I've been sat infront of my computer for almost a day, and still can't see
why this wont work. I've made a raster that has split the world into 49
different regions. I'm then going to input the suitabilit

Re: raster package options maxmemory and chunksize

Tim Howard — 2012-05-17T12:05:32

Some updates for the list archives:
 
Using  
setOptions(maxmemory = 1e+09, chunksize=1e+08)
reduced predict time to about 10 hours 50 min on my system. Good!  I also spent the night converting all the GRIDS to native raster format and am trying the next run with that setup to see if read times differ significantly. 
 
I'll report back, and then also see if I can increase chunksize even further. 
 
The raster package truly simplifies things  -- I love it!
 
R version 2.15.0 (2012-03-30)
Platform: x86_64-pc-mingw32/x64 (64-bit)
 
locale:
[1] LC_COLLATE=English_United States.1252 
[2] LC_CTYPE=English_United States.1252   
[3] LC_MONETARY=English_United States.1252
[4] LC_NUMERIC=C                          
[5] LC_TIME=English_United States.1252    
 
attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     
 
other attached packages:
[1] rgdal_0.7-8        randomForest_4.6-6 raster_1.9-82      sp_0.9-98         
 
loaded via a namespace (and not attached):
[1] grid_2.15

Re: Adding additional items to plot

Oscar Perpiñán Lamigueiro — 2012-05-17T09:04:26


Hello,

Try this approach using rasterVis and latticeExtra:

library(raster)
library(rasterVis)
library(maptools) ##for readShapeLines
library(latticeExtra)

## Get the data from gadm and data-vis
download.file('http://www.gadm.org/data/shp/ETH_adm.zip', '~/temp/ETH_adm.zip')
unzip('~/temp/ETH_adm.zip', exdir='~/temp')
download.file('http://www.diva-gis.org/data/msk_alt/ETH_msk_alt.zip', '~/temp/ETH_msk_alt.zip')
unzip('~/temp/ETH_msk_alt.zip', exdir='~/temp')
download.file('http://www.diva-gis.org/data/wat/ETH_wat.zip', '~/temp/ETH_wat.zip')
unzip('~/temp/ETH_wat.zip', exdir='~/temp')

## Read the files
proj <- CRS(' +proj=longlat +ellps=WGS84')
ethAdm <- readShapeLines('~/temp/ETH_adm1.shp', proj4string=proj)
ethWat <- readShapeLines('~/temp/ETH_water_lines_dcw.shp', proj4string=proj)
ethAlt <- raster('~/temp/ETH_msk_alt')

## and plot using the "+.trellis" and "layer" functions
## from latticeExtra
levelplot(ethAlt, par.settings=GrTheme) +
  layer(sp.lines(ethWat, col='blue', lwd=0.6)) +
  layer(sp.line

Coordinate transformation for directions?

Barry Rowlingson — 2012-05-17T08:55:03

spTransform does a great job of transforming points from one CRS to
another, but what's the best way of transforming directions?

For example, the wind at a point is blowing from the SW in lat-long
coordinates. That would be -3pi/4 radians if you consider lat-long as
a cartesian coordinate system. But suppose you want to work in
Ordnance Survey grid references. These coords are rotated slightly
compared to epsg:4326, and the rotation is different at different
points. Lat-long North is never exactly straight up on OSGrid maps.

 My thinking is to consider a small vector (how small?) from (x,y) to
(x+dx,y+dy) in the required direction in lat-long, then transform
those two points to OSGrid, then compute the angle with the atan2
function.

 Apart from wondering how small the vector needs to be, have I missed
any other subtleties? Has anyone written this already?

 Note these are point-wise directions, not bearings of A from B. Great
circles are not involved!

Barry

Re: rgdal package for mac

Gabriele Cozzi — 2012-05-17T08:33:21

Hi Roger,
sorry for troubling and thanks for the very quick answer.

I did follow the guidelines from the CRAN page
# setRepositories(ind=1:2)
# install.packages('rgdal'),
but it did not work.

I was successful, however (with the same to lines of code), after 
shutting down and restating R...

============================================================

Gabriele Cozzi

Institute of Evolutionary Biology and Environmental Studies
Zurich University
Winterthurerstr. 190
8057 Zurich - Switzerland
E-mail: gabriele.cozzi< at >uzh.ch
Phone: ++41(0)44 6356116
Fax: ++41(0)44 6355711
http://www.ieu.uzh.ch
http://african-carnivores.com

Botswana Predator Conservation Trust
Private Bag 13
Maun - Botswana
E-mail: gab.cozzi< at >gmail.com
Mobile: +26774259312
http://www.bpctrust.org


On 17.05.2012 10:06, Roger Bivand wrote:

Re: rgdal package for mac

Roger Bivand — 2012-05-17T08:06:20


Did you read the CRAN package page for rgdal? If you did, did you do what 
it said:

setRepositories(ind=1:2)
install.packages('rgdal')

to install from CRAN Extras?

Don't jump to conclusions, do try to use the resources available to reach 
well grounded solutions.

We still hope that rgdal and rgeos will be available as regular CRAN 
packages for OSX, but this seems to have stranded (again).

Roger

rgdal package for mac

Gabriele Cozzi — 2012-05-17T08:05:47

Dear all,

I am having issues installing the package "rgdal" on a Mac OS X (Lion).

# install.packages("rgdal")

actually says that there is no available package for R 2.15.0 (the 
version I am using at the moment). Yet I can successfully load the rgdal 
package on my old PC (also running R 2.15.0). Does it means rgdal is not 
implemented for mac? I find it hard to believe.

Any help would be appreciated.

Best,
Gabriele

============================================================

Gabriele Cozzi

Institute of Evolutionary Biology and Environmental Studies
Zurich University
Winterthurerstr. 190
8057 Zurich - Switzerland
E-mail: gabriele.cozzi< at >uzh.ch
Phone: ++41(0)44 6356116
Fax: ++41(0)44 6355711
http://www.ieu.uzh.ch
http://african-carnivores.com

Botswana Predator Conservation Trust
Private Bag 13
Maun - Botswana
E-mail: gab.cozzi< at >gmail.com
Mobile: +26774259312
http://www.bpctrust.org


On 17.05.2012 09:15, Roger Bivand wrote:

Re: trouble with spTransform

Roger Bivand — 2012-05-17T07:15:26


Your points are 750m east and 600m north of the given projection origin. 
If you meant km, set the units by adding +units=km to the proj4string. 
However, this gives me:

          coordinates Loc
1  (101.24, 5.44994)   a
2 (101.295, 5.56496)   b

which is not on Phuket, so something else is wrong here. Do you know 
definitely that the input coordinates are UTM zone 47? Reversing the 
eastings and northings gets a bit closer:

          coordinates Loc
1  (99.9308, 6.7677)   a
2 (100.046, 6.82295)   b

but isn't there yet (reversed coordinate order does happen). What was the 
source of the coordinates?

Roger

trouble with spTransform

Hodgess, Erin — 2012-05-17T06:02:10

Hello again.

I'm having a little bit of trouble with spTransform (probably I'm doing it wrong), but here is what I'm doing:


Original data frame with UTM locations near Phuket, Thailand:
  Loc    east   north
1   a 748.168 602.861
2   b 754.302 615.609

Set up a spatial data frame
         coordinates Loc
1 (748.168, 602.861)   a
2 (754.302, 615.609)   b
           coordinates Loc
1 (94.518, 0.00543748)   a
2 (94.518, 0.00555246)   b

Now the answers should be about 8 degrees north and 98 degrees East.

Any suggestions would be most appreciated.

Sincerely,
Erin


Erin M. Hodgess, PhD
Associate Professor
Department of Computer and Mathematical Sciences
University of Houston - Downtown
mailto: hodgesse< at >uhd.edu


[[alternative HTML version deleted]]

Re: converting from UTM to Decimal degrees

Hodgess, Erin — 2012-05-17T05:18:20

thanks!

I knew about spTransform but not project.

Sincerely,
Erin


Erin M. Hodgess, PhD
Associate Professor
Department of Computer and Mathematical Sciences
University of Houston - Downtown
mailto: hodgesse< at >uhd.edu



-----Original Message-----
From: Michael Sumner [mailto:mdsumner< at >gmail.com]
Sent: Thu 5/17/2012 12:15 AM
To: Hodgess, Erin
Cc: r-sig-geo< at >r-project.org
Subject: Re: [R-sig-Geo] converting from UTM to Decimal degrees
 
See functions ?project and ?spTransform in the rgdal package. The
former is simpler but is only for transforming  a matrix of
coordinates from or to longitude / latitude, and no datum
transformation is possible. The latter can transform Spatial* objects
arbitrarily from or or to any PROJ.4 coordinate system and includes
datum transformations.


On Thu, May 17, 2012 at 2:47 PM, Hodgess, Erin <HodgessE< at >uhd.edu> wrote:

Re: converting from UTM to Decimal degrees

Michael Sumner — 2012-05-17T05:15:47

See functions ?project and ?spTransform in the rgdal package. The
former is simpler but is only for transforming  a matrix of
coordinates from or to longitude / latitude, and no datum
transformation is possible. The latter can transform Spatial* objects
arbitrarily from or or to any PROJ.4 coordinate system and includes
datum transformations.


On Thu, May 17, 2012 at 2:47 PM, Hodgess, Erin <HodgessE< at >uhd.edu> wrote:

converting from UTM to Decimal degrees

Hodgess, Erin — 2012-05-17T04:47:10

Dear R Sig Geo People:

Is there a function to convert UTM coordinates to Decimal degrees, please?  This is just for a few points; nothing big.

I know that this is a relatively easy one, but thought I'd check before I wrote up something.

Thanks,
Erin



Erin M. Hodgess, PhD
Associate Professor
Department of Computer and Mathematical Sciences
University of Houston - Downtown
mailto: hodgesse< at >uhd.edu


[[alternative HTML version deleted]]

Learning patterns from temporal / spatial data

2012-05-17T04:32:45

Hi,

  Iam a newbie to R as well as analysis of temporal / spatial data.



I have a set of GPS coordinates each associated with time,  speed as well as vehicle id.



Sample data points are as below  :
North_Cod

East_Cod

Date1

Time1

Dt1Tm1

Veh_Id

Veh_Speed

77.542016

28.437566

3/12/2012 0:00

0

3/12/2012 0:00

482

0

77.542016

28.437566

3/12/2012 0:00

0

3/12/2012 0:00

482

0

77.542016

28.437566

3/12/2012 0:00

0

3/12/2012 0:00

482

0

77.585216

28.42415

3/12/2012 0:00

0

3/12/2012 0:00

434

0

77.585216

28.42415

3/12/2012 0:00

0

3/12/2012 0:00

434

0








I want to do following

1. Visual analysis of this data

2. Learn some patterns on this data by clustering



a)  Which is the best package / packages in R for above

b)  Are there any ready examples which I can use readily for above



A high level approach on how I could go about for above would be very helpful.





Thanks for all help,

Regards,

Arthi







Please do not print this email unless it is absolutely necessa

Adding additional items to plot

Adam Sparks — 2012-05-17T02:54:56

Dear all,
I'm attempting to produce some finished product maps using R and have been
struggling with an issue.

I want to use a hillshade object created from a DEM using the hillShade()
function as a backdrop for my raster data that I'm presenting.

Behind that I want to use lightsteelblue2 as the color where there is water.

On top of this I want to plot the country outlines for the areas of
interest.

So in summary I have:
 the first layer, blue for water
 the second layer, hillshade DEM
 the third layer, my raster data I'm presenting
 the fourth layer, the country outline

Two problems arise. The first is, if I create an empty plot using plot(type
= 'n') I don't have the correct margins for a legend.

So I can just plot my raster data once, then put the blue rectangle on top
of that and I have the proper margins.  Not elegant but it works.

The second issue that I've not been able to sort out is that when I plot
the country outline it doesn't respect the plot area. It acts as if there
is no margin there.

Re: Calculating area of region with a raster

Anthony Fischbach — 2012-05-16T22:56:43

Rhys,
using the NASA special sensor microwave imagery sea ice data set with its
constant area pixel 
(PixelArea = 25x25 ## in km) geometry I was able to solve this problem by
creating a raster with pixel values = 1 for pixels with ice and NA for
pixels without ice then I summed the values and multiplied by the PixelArea.
sum(values(raster.Ice), na.rm=T)* PixelArea
hth-

-----
Tony Fischbach, Wildlife Biologist
Walrus Research Program
Alaska Science Center
U.S. Geological Survey
4210 University Drive
Anchorage, AK 99508-4650

AFischbach< at >usgs.gov
http://alaska.usgs.gov/science/biology/walrus
--
View this message in context: http://r-sig-geo.2731867.n2.nabble.com/Calculating-area-of-region-with-a-raster-tp7561970p7562866.html
Sent from the R-sig-geo mailing list archive at Nabble.com.

raster package options maxmemory and chunksize

Tim Howard — 2012-05-16T16:58:12

all - 
 
I'm trying to tweak my setup of the raster package to increase the speed of a predict on a randomForest object. 
 
My system: Windows 7 64-bit with 12GB RAM using R x64 2.15.0, raster version 1.9-82
 
I'm running the predict on a rasterStack of about 45 layers. 
 
Question #1. Does chunksize refer to the chunk of a single layer in the stack or the entire stack? To rephrase, do I have to be sure there is room to take a chunk consisting of chunksize from every layer and thus need room for chunksize * 45?
 
The default is 1e+07, which I read as 10MB.  Is that right?
 
Question #2. Is there a rule of thumb for how much buffer to leave between chunksize and maxmemory? The default is 1+e08, which I read as 100MB.
 
I just tried 
 
setOptions(maxmemory = 4e+09, chunksize=1e+09)
 
and the predict never got past 0% completion before bailing with a memory error. 
 
Dropping it down to
 
setOptions(maxmemory = 1e+09, chunksize=1e+08)
 
is looking like it has promise, with RGui using about 1-2GB. 
 
Using the d

Re: reading a large TXT file and save in chunk

Robert J. Hijmans — 2012-05-15T22:09:58

Hi Gianni,

Here is a link:  http://comments.gmane.org/gmane.comp.lang.r.general/255857

In my example, replace  "i <- i + length(d)"  with "write.table( )"

Robert



On Tue, May 15, 2012 at 1:15 PM, Sean O'Riordain <seanpor< at >acm.org> wrote:


[[alternative HTML version deleted]]