Subscribe to R-sig-geo feed
This is an archive for R-sig-geo, not a forum. It is not possible to post through Nabble - you may not start a new thread nor follow up an existing thread. If you wish to post, but are not subscribed to the list through the list homepage, subscribe first (through the list homepage, not via Nabble) and post from your subscribed email address. Until 2015-06-20, subscribers could post through Nabble, but policy has been changed as too many non-subscribers misunderstood the interface.
Updated: 35 min 11 sec ago

Paired (before-after comparison) t-test with spatial locations

Fri, 06/23/2017 - 20:07
Hello R-sig-geo list,

We have soil chemistry data before and after fire at specific spatial
locations.

Is there a "paired t-test" equivalent that can incorporate geospatial
coordinates and hence adjust for spatial dependence in the data?

Apologies in advance if I have missed the obvious in *spdep* and *spatstat*,
but I am new to *R* and would appreciate your advice.

Dora Pearce
Honorary: Faculty of Science and Technology
Federation University Australia
University Dr, Mount Helen VIC 3350
https://federation.edu.au/faculties-and-schools/faculty-of-science-and-technology

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: Problems With Cartography 2 - Help

Fri, 06/23/2017 - 05:44
Thank you Chris,

I found a way to connect the three islands in a previous R-sig-Geo answer
from Roger Bivand (
https://stat.ethz.ch/pipermail/r-sig-geo/2013-July/018868.html). I followed
the code given in this answer and it seems to work fine (I plotted it and
it looks ok).

x.nb <- poly2nb(carto2012u2)
summary(x.nb)

### 1) fIND AREAS WITH NO NEIGHBORS
no_neighs <- which(card(x.nb) == 0) # returns set cardinality.

### To check that you get the same result using counts of graph components:

nc <- n.comp.nb(x.nb)
table(nc$comp.id)
no_neighs_comps <- which(unname(table(nc$comp.id)==1))
match(no_neighs_comps, nc$comp.id)

### 2) For each island find nearest neighbour

k1nb <- knn2nb(knearneigh(coordinates(carto2012u2), k=1,
                          longlat=!is.projected(carto2012u2)))

### 3) Assign a neighbour relationship between the island and its nearest
neighbour

is.symmetric.nb(x.nb, force=TRUE)
nb1 <- x.nb

res <- k1nb[no_neighs]
nb1[no_neighs] <- res
attr(nb1, "sym") <- is.symmetric.nb(nb1, force=TRUE)

# re-assign the now incorrect symmetry attribute

nb1

nb2 <- make.sym.nb(nb1)
is.symmetric.nb(nb2, force=TRUE)
nb2

table(n.comp.nb(nb2)$comp.id)

###

veins_inla <- nb2INLA("carto2012_nb.inla", nb2)

plot (carto2012u2, axes=F)
coords <- coordinates(carto2012u2)
plot(nb2, coords, add=T, col="red")

Hope this would be helpful for others with the same issue.


Andrés


On Thu, May 11, 2017 at 3:25 PM, chris english <
[hidden email]> wrote:

> Andres,
>
> I'm thinking you want to 'unsimplify' the topology prior to your poly2nb
> by a slight negative gBuffer on the Galapagos polygons,
> reduce the Galapagos polygons by the buffer and un-share the boundaries.
> I'm trying some code, but looking at a mix of these as
> guidance:
>
> https://www.rdocumentation.org/packages/rgeos/versions/0.
> 3-22/topics/gBuffer
> https://gist.github.com/mstrimas/1b4a4b93a9d4a158bce4
> https://gis.stackexchange.com/questions/93096/how-to-
> perform-a-true-gis-clip-of-polygons-layer-using-a-polygon-layer-in-r
>
> Saludos,
> Chris
>
> On Wed, May 10, 2017 at 4:27 AM, Andrés Peralta <[hidden email]>
> wrote:
>
>> Hi to everyone,
>>
>> I`m working with the second administrative level cartography from Ecuador
>> (available in: http://www.ecuadorencifras.gob.ec//documentos/web-inec/
>> Cartografia/2015/Clasificador_Geografico/2012/SHP.zip
>> <http://www.ecuadorencifras.gob.ec//documentos/web-inec/Cartografia/2015/Clasificador_Geografico/2012/SHP.zip>).
>> The shape file is
>> called nxcantones.shp. I`ve been having a lot of problems opening and
>> working with this cartography in R; but i can open it in Q-GIS and ARCGIS
>> without any problem (I have opened it in both programs and saved it
>> again).
>> As the cartography has various objects with the same ID - aparently the
>> same areas but repeated; we had to run the following sintax in order to
>> have one ID in each area:
>>
>> *carto2012 <- readOGR("CANTONES2012/2012CLEAN.shp", "2012CLEAN",
>> stringsAsFactors=F) *
>> *proj4string(carto2012) <- CRS("+proj=utm +zone=17 +south +datum=WGS84
>> +units=m +no_defs")*
>>
>> *carto2012x <- unionSpatialPolygons(carto2012, IDs=carto2012$DPA_CANTON)*
>> *proj4string(carto2012x) <- CRS("+proj=utm +zone=17 +south +datum=WGS84
>> +units=m +no_defs")*
>>
>>
>> *datos <- data.frame(ID=row.names(carto2012x), stringsAsFactors = F)*
>> *row.names(datos) <- row.names(carto2012x)*
>> *carto2012x2 <- SpatialPolygonsDataFrame(carto2012x, data=datos)*
>>
>> *carto2012 <- carto2012x2 *
>>
>> For our work, we had to make a neighborhood matrix of the cartography.
>>
>> *x.nb <- poly2nb(carto2012u2) # U2 is the cartography after the recode of
>> some areas that had to be joined.*
>>
>> *summary(x.nb)*
>>
>> This seems to work fine. The problem is we need to connect the three
>> island
>> areas (the Galapagos) but when we try to do so, it connects many areas
>> along all the cartography. I've tried to do it manually using *edit.nb*
>> but
>> it does not work. I´ve also tried the following sintaxes:
>>
>> x.nb1 <- x.nb
>> which(card(x.nb1)==0) #to discover the id of these areas without
>> connections
>> id <- function(x){which(carto2012u2$ID==x)}
>>
>> x.nb1[[id(2001)]] = unique(as.integer(sort(c(x.nb1[[id(2001)]],
>> id(2002)))))
>> x.nb1[[id(2002)]] = unique(as.integer(sort(c(x.nb1[[id(2002)]],
>> id(2001)))))
>>
>> x.nb1[[id(2001)]] = unique(as.integer(sort(c(x.nb1[[id(2001)]],
>> id(2003)))))
>> x.nb1[[id(2003)]] = unique(as.integer(sort(c(x.nb1[[id(2003)]],
>> id(2001)))))
>>
>> x.nb1[[id(2002)]] = unique(as.integer(sort(c(x.nb1[[id(2002)]],
>> id(2003)))))
>> x.nb1[[id(2003)]] = unique(as.integer(sort(c(x.nb1[[id(2003)]],
>> id(2002)))))
>>
>> ####
>> x.nb1[[id("2001")]] = unique(as.integer(sort(c(x.nb1[[id("2001")]],
>> id("2002")))))
>> x.nb1[[id("2002")]] = unique(as.integer(sort( id("2001"))))
>>
>> x.nb1[[id("2001")]] = unique(as.integer(sort(c(x.nb1[[id("2001")]],
>> id("2003")))))
>> x.nb1[[id("2003")]] = unique(as.integer(sort( id("2001"))))
>>
>> x.nb1[[id("2002")]] = unique(as.integer(sort(c(x.nb1[[id("2002")]],
>> id("2003")))))
>> x.nb1[[id("2003")]] = unique(as.integer(sort( id("2002"))))
>>
>> Any ideas or suggestions? Your help will really be appreciated.
>>
>>
>> --
>>
>> *             Andrés Peralta*
>>         Pre-doctoral Researcher
>>      GREDS/EMCONET / ASPB
>> <https://www.upf.edu/greds-emconet/en/>   <http://aspb.cat/>
>>
>>         [[alternative HTML version deleted]]
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> [hidden email]
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
>
>

--

*             Andrés Peralta*
        Pre-doctoral Researcher
     GREDS/EMCONET / ASPB
<https://www.upf.edu/greds-emconet/en/>   <http://aspb.cat/>

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Handling Spatial Autocorrelation and Multicolinearity when modelling altitude of forestline

Fri, 06/23/2017 - 04:46

Dear all,

 

Earlier I asked this question on GRASS ML [1], but maybe it is more suited here:

Zofie and me are trying to mode altitude of forest line in Fensocandia using GRASS and R.

 

We have a couple of 100k pixels that we assume to represent forest line and now want to model their altitude with explanatory variables (latitude, terrain, temperature, precipitation and the like). In a next step we want to use projected data (climate scenarios) in the model in order to predict possible effects of climate change. But now we are a bit unsure about what modeling technique to use.

 

Our data to be modeled is zero-inflated (from 0 - ~1150m) with a significant amount of spatial autocorrelation. And also the explanatory variables are spatially auto-correlated and have a lot of collinearity.

We have been looking at (amongst others):

  • Regression kriging, but we have doubts that R will be able to handle the amount of data (even on a highmem server), pluss that we are unsure if we can replace current with future climate data in such a model
  • GLS, but also here we face problems with excessive resource consumption in R for e.g. accounting for spatial autocorrelation and on a subsample residuals still showed spatial autocorrelation.

 

Can anyone recommend other types of models or R packages we should look at, that can handle spatial autocorrelation and multicolinearity?

We would be glad for any hint also on where to look for more information (articles, textbooks…)?

 

Thanks for helping in advance!

 

Cheers

Stefan

 

1: https://www.mail-archive.com/[hidden email]/msg34297.html


_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: raster - unrotate?

Thu, 06/22/2017 - 18:26
On Thu, 22 Jun 2017 at 21:28 Ben Tupper <[hidden email]> wrote:

> Hi,
>
> Wow!  A treasure from the past!
>
> Hmmm.   I see what you mean; it might be best to snip-clip-and-smush from
> the native presentation.  We are testing out the performance of the dineof
> function in the sinkr package (
> https://github.com/marchtaylor/sinkr/blob/master/R/dineof.R) to
> interpolate patchy chlorophyl data.
>
> Thanks for the tip!
>
>
For what it's worth, I'm happy to help, there's no one size fits all here.
The dateline is one of those "easy problems" that does not have a general
fix and will bite in many different contexts. :)

The best distillation of my tricks is in this project, but it does depend
on your workflow and the actual data in use.

https://github.com/hypertidy/tabularaster

If you have performance issues raster's cell abstraction  will help, but
they are a bit old-school when it comes to multi-part geometries. (It's
identical to standard database indexing concepts as are all "spatial"
optimization tricks)

(I see a desperate need for a general API for this kind of problem so that
we can build tools from a general framework, which I'm working on - that's
some kind of excuse for why this and related projects are quite raw and
unfinished.)

Cheers, Mike.

Cheers,
> Ben
>
> On Jun 22, 2017, at 4:46 AM, Michael Sumner <[hidden email]> wrote:
>
>
> It used to do the inverse, and I preserved a copy of the old behaviour
> here:
>
>
> https://github.com/AustralianAntarcticDivision/raadtools/blob/master/R/utils.R#L9
>
> You're probably best to learn from how it approaches it rather than just
> use it, since it's not a general capability, just works for those very
> specific cases.
>
> I just tested it and it still seems to work fine, I used
>
>
> oceandata.sci.gsfc.nasa.gov/MODISA/Mapped/Monthly/9km/chlor/A20021822002212.L3m_MO_CHL_chlor_a_9km.nc
>
>
> obtainable with
>
>
> https://oceandata.sci.gsfc.nasa.gov/cgi/getfile/A20021822002212.L3m_MO_CHL_chlor_a_9km.nc
>
>
> If you're extracting points you might as well just convert them into the
> "Atlantic" range, but it's painful to get that right for polygons or lines
> (and very hard to generalize once you get it working anyway).
>
> For simple regions like the one you show I'd probably split it into two
> extents() and use those - but depends on your workflow, all these options
> are useful here and there.
>
> Cheers, Mike.
>
>
> On Thu, 22 Jun 2017 at 18:30 Ben Tupper <[hidden email]> wrote:
>
>> Hello,
>>
>> We have rasters that span [-180, 180] from NASA's Ocean Color (
>> https://oceancolor.gsfc.nasa.gov/)  datasets.  We are trying to extract
>> a region that spans 100E to 90W, that is 100 to -90.  The region 'wraps'
>> across the edges as shown by the plot at the address below.
>>
>> https://dl.dropboxusercontent.com/u/8433654/sst.png
>>
>>
>> library(raster)
>> uri <- '
>> https://oceandata.sci.gsfc.nasa.gov:443/opendap/MODISA/L3SMI/2016/001/A20160012016032.L3m_R32_SST_sst_9km.nc
>> <https://oceandata.sci.gsfc.nasa.gov/opendap/MODISA/L3SMI/2016/001/A20160012016032.L3m_R32_SST_sst_9km.nc>
>> '
>> R <- raster::raster(uri, varname = 'sst')
>>
>> R
>> #class       : RasterLayer
>> #dimensions  : 180, 360, 64800  (nrow, ncol, ncell)
>> #resolution  : 1, 1  (x, y)
>> #extent      : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
>> #coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
>> #data source : in memory
>> #names       : layer
>> #values      : 1.482572e-05, 0.9999928  (min, max)
>>
>> plot(R)
>> lines(c(180, 100, 100, 180), c(80,80,0,0))
>> lines(c(-180,-90,-90,-180), c(80,80,0,0))
>>
>> I see that there is the nice raster::rotate() function to rotate raster
>> coordinates from [0,360] to [-180,180]. That would make extracting the
>> region super easy if the inverse were available.  Is there an equivalent
>> way to implement the inverse or raster::rotate()?  I could dig into the
>> source of raster::rotate() to see how to code one up, but I hate like heck
>> to reinvent the wheel.
>>
>> Thanks!
>> Ben
>>
>>
>> Ben Tupper
>> Bigelow Laboratory for Ocean Sciences
>> 60 Bigelow Drive, P.O. Box 380
>> East Boothbay, Maine 04544
>> http://www.bigelow.org
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> [hidden email]
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
> --
> Dr. Michael Sumner
> Software and Database Engineer
> Australian Antarctic Division
> 203 Channel Highway
> Kingston Tasmania 7050 Australia
>
>
> Ben Tupper
> Bigelow Laboratory for Ocean Sciences
> 60 Bigelow Drive, P.O. Box 380
> East Boothbay, Maine 04544
> http://www.bigelow.org
>
>
>
> -- Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: sf 0.5-0 binary for macOS SIerra?

Thu, 06/22/2017 - 08:11
On Thu, Jun 22, 2017 at 8:58 AM, Roger Bivand <[hidden email]> wrote:

> On Thu, 22 Jun 2017, Roger Bivand wrote:
>
> See this thread (not only Roy's brain that kicks in more slowly ... ):
>
> https://stat.ethz.ch/pipermail/r-sig-mac/2017-June/012429.html


Thank you, I will try that and report back.

Kent

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: Thresholds & reclassify raster

Thu, 06/22/2017 - 07:59
Many thanks Ben Tupper!
Best Rgds John

On Thu, Jun 22, 2017 at 1:47 PM, Ben Tupper <[hidden email]> wrote:

> Hi,
>
> Your attachment didn't make it through.
>
> I'm not sure about your first option, but your second can use a
> combination of quantile() and findInterval() to create a classified
> raster.  Perhaps like this...
>
>
> library(raster)
> r <- raster()
> r <- setValues(r, runif(ncell(r), 40, 90))
>
> pp <- quantile(r, c(0, 0.15, 0.85))
> # pp
> #       0%      15%      85%
> # 40.00199 47.64569 82.50751
>
> ix <- findInterval(getValues(r), pp)
>
> classified_r <- setValues(r, ix)
>
> # > r
> # class       : RasterLayer
> # dimensions  : 180, 360, 64800  (nrow, ncol, ncell)
> # resolution  : 1, 1  (x, y)
> # extent      : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
> # coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
> # data source : in memory
> # names       : layer
> # values      : 40.00199, 89.99892  (min, max)
>
> # > classified_r
> # class       : RasterLayer
> # dimensions  : 180, 360, 64800  (nrow, ncol, ncell)
> # resolution  : 1, 1  (x, y)
> # extent      : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
> # coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
> # data source : in memory
> # names       : layer
> # values      : 1, 3  (min, max)
>
>
> Cheers,
> Ben
> > On Jun 22, 2017, at 4:04 AM, John Wasige <[hidden email]> wrote:
> >
> > ​Dear all​
> > ,
> > I am requesting for your help to reclassify the gpp.tif (attached) raster
> > in R. Not sure my script here is doing the right thing.
> >
> > ​#Data structure​
> >
> > class       : RasterBrick
> > dimensions  : 1243, 1409, 1751387, 1  (nrow, ncol, ncell, nlayers)
> > resolution  : 0.00225804, 0.00225804  (x, y)
> > extent      : -25.53282, -22.35125, 14.54232, 17.34906  (xmin, xmax,
> ymin, ymax)
> > coord. ref. : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84
> +towgs84=0,0,0
> > data source : L:\NDVI_WGS84\gpp.tif
> > names       :        gpp
> > min values  : -0.2176632
> > max values  :   2.989293
> >
> >
> > ​###​
> >
> > ​I want to ​
> > use the statistical distribution of gpp.tif values set threshold for
> change
> > assignment. The thresholds should be set by the following two
> > ​ ​
> > options:
> >
> > (Option1) => use  mean and StD and assign all pixel values  < mean - StD
> as
> > negative change (class 1) and all values  > mean + StD as positive
> > change (class
> > 3); all values in between would be assigned 'neutral' for the final
> > tabulation (class2).
> >
> > (option2) => what may come to almost the same result is to do the
> > cumulative histogram of the calculated GPPproxy pixels and assign all
> > pixels below the 15 percentile to negative change (class1) and all pixel
> > above the 85 percentile to positive change (class3) and again the range
> in
> > between (15 to 85 percentile) would be assigned neutral (class2).
> >
> > #R script
> >
> > # classification
> > #Option1_mean_sd
> > gpp <- brick('L:/NDVI_WGS84/gpp.tif')
> > hist(as.matrix(GPP))
> > plot(GPPproxy)
> > m <- cellStats(gpp, 'mean')
> > sd <- cellStats(gpp, 'sd')
> > v <- m-sd
> > gppnew <- gpp-v
> > plot(gppnew)
> > gppnew[gppnew>0] <- 3
> > gppnew[gppnew<0] <- 1
> > gppnew[gppnew==0] <- 2
> > plot(gppnew)
> > writeRaster(gppnew, filename=paste("L:/NDVI_WGS84/
> gppnew_mean_sd",sep=''),
> > format="GTiff", overwrite=TRUE)
> >
> > #Option2_percentile
> > max<- cellStats(gpp, 'max')
> > gppnew2 <-(gpp/max)*100
> > min <- cellStats(gpp, 'min')
> > m <- matrix(c(min, 15, 1,
> >              15, 85, 2,
> >              85, Inf, 3), ncol=3, byrow=TRUE)
> > gppnew_percentile <- reclassify(gppnew2, m)
> >
> > hist(as.matrix(gppnew_percentile))
> > writeRaster(gppnew_percentile,
> > filename=paste("L:/NDVI_WGS84/gppnew_percentile",sep=''),
> > format="GTiff", overwrite=TRUE)
> > plot(gppnew_percentile)
> >
> > ​Thanks for your help
> >
> >       [[alternative HTML version deleted]]
> >
> > _______________________________________________
> > R-sig-Geo mailing list
> > [hidden email]
> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
> Ben Tupper
> Bigelow Laboratory for Ocean Sciences
> 60 Bigelow Drive, P.O. Box 380
> East Boothbay, Maine 04544
> http://www.bigelow.org
>
>
>
>

--
John Wasige
"There are no REGRATES in LIFE, just lessons (Jennifer Aniston)”

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: sf 0.5-0 binary for macOS SIerra?

Thu, 06/22/2017 - 07:58
On Thu, 22 Jun 2017, Roger Bivand wrote:

> Please also follow:
>
> https: //github.com/edzer/sfr/issues/393
>
> and consider trying:
>
> http://www.kyngchaos.com/software/frameworks

See this thread (not only Roy's brain that kicks in more slowly ... ):

https://stat.ethz.ch/pipermail/r-sig-mac/2017-June/012429.html

Roger

>
> for a source sf install until things stabilise.
>
> Roger
>
> On Thu, 22 Jun 2017, Kent Johnson wrote:
>
>>  Is there a simple features 0.5-0 binary for macOS Sierra (10.12)? I am
>>  having trouble installing it from source. CRAN shows
>>  OS X El Capitan binaries:r-release: sf_0.4-3.tgz
>>  OS X Mavericks binaries:r-oldrel: sf_0.5-0.tgz
>>
>>  Can I use the Mavericks binary?
>>
>>  Thanks,
>>  Kent
>>
>>   [[alternative HTML version deleted]]
>>
>>  _______________________________________________
>>  R-sig-Geo mailing list
>>  [hidden email]
>> https: //stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>
>
--
Roger Bivand
Department of Economics, Norwegian School of Economics,
Helleveien 30, N-5045 Bergen, Norway.
voice: +47 55 95 93 55; e-mail: [hidden email]
Editor-in-Chief of The R Journal, https://journal.r-project.org/index.html
http://orcid.org/0000-0003-2392-6140
https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Roger Bivand
Department of Economics
Norwegian School of Economics
Helleveien 30
N-5045 Bergen, Norway

Re: Problems with using dismo package on HPC

Thu, 06/22/2017 - 07:08
It might be because the file is absent.
You need to trace through from the top, what file path is passed to raster,
then try reading that directly with raster() and ensure the file exists.
It's possibly a raster object with a link to a file that no longer exists.
Raster is a complex polymorphic beast, and can be a thin brittle shell held
together by a file name, with no safety net.

If you don't have the original code it's harder, but also hard to help
without access to your system.

Cheers, Mike

On Wed, 21 Jun 2017, 07:30 Wall, Wade A ERDC-RDE-CERL-IL CIV, <
[hidden email]> wrote:

> Hi all,
>
> I am trying to use the dismo package on an HPC and am running into some
> problems. The code that I used worked fine until a few days ago, then I
> started getting the error message
>
> Error in rgdal::getRasterData(con, offset = offs, region.dim = c(1, nc),  :
> Failure during raster IO
>
> Not sure how to proceed on this. The code works on my desktop computer,
> and like I said, worked until a few days ago on the HPC. I have read/write
> privileges to the outfolder location,  all the environmental layers are the
> same extent and resolution, and the spatial points are in the same
> coordinate system.
>
> Any help would be greatly appreciated.
> Wade
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> --
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: Thresholds & reclassify raster

Thu, 06/22/2017 - 06:47
Hi,

Your attachment didn't make it through.

I'm not sure about your first option, but your second can use a combination of quantile() and findInterval() to create a classified raster.  Perhaps like this...


library(raster)
r <- raster()
r <- setValues(r, runif(ncell(r), 40, 90))

pp <- quantile(r, c(0, 0.15, 0.85))
# pp
#       0%      15%      85%
# 40.00199 47.64569 82.50751

ix <- findInterval(getValues(r), pp)

classified_r <- setValues(r, ix)

# > r
# class       : RasterLayer
# dimensions  : 180, 360, 64800  (nrow, ncol, ncell)
# resolution  : 1, 1  (x, y)
# extent      : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
# coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
# data source : in memory
# names       : layer
# values      : 40.00199, 89.99892  (min, max)

# > classified_r
# class       : RasterLayer
# dimensions  : 180, 360, 64800  (nrow, ncol, ncell)
# resolution  : 1, 1  (x, y)
# extent      : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
# coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
# data source : in memory
# names       : layer
# values      : 1, 3  (min, max)


Cheers,
Ben
> On Jun 22, 2017, at 4:04 AM, John Wasige <[hidden email]> wrote:
>
> ​Dear all​
> ,
> I am requesting for your help to reclassify the gpp.tif (attached) raster
> in R. Not sure my script here is doing the right thing.
>
> ​#Data structure​
>
> class       : RasterBrick
> dimensions  : 1243, 1409, 1751387, 1  (nrow, ncol, ncell, nlayers)
> resolution  : 0.00225804, 0.00225804  (x, y)
> extent      : -25.53282, -22.35125, 14.54232, 17.34906  (xmin, xmax, ymin, ymax)
> coord. ref. : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
> data source : L:\NDVI_WGS84\gpp.tif
> names       :        gpp
> min values  : -0.2176632
> max values  :   2.989293
>
>
> ​###​
>
> ​I want to ​
> use the statistical distribution of gpp.tif values set threshold for change
> assignment. The thresholds should be set by the following two
> ​ ​
> options:
>
> (Option1) => use  mean and StD and assign all pixel values  < mean - StD as
> negative change (class 1) and all values  > mean + StD as positive
> change (class
> 3); all values in between would be assigned 'neutral' for the final
> tabulation (class2).
>
> (option2) => what may come to almost the same result is to do the
> cumulative histogram of the calculated GPPproxy pixels and assign all
> pixels below the 15 percentile to negative change (class1) and all pixel
> above the 85 percentile to positive change (class3) and again the range in
> between (15 to 85 percentile) would be assigned neutral (class2).
>
> #R script
>
> # classification
> #Option1_mean_sd
> gpp <- brick('L:/NDVI_WGS84/gpp.tif')
> hist(as.matrix(GPP))
> plot(GPPproxy)
> m <- cellStats(gpp, 'mean')
> sd <- cellStats(gpp, 'sd')
> v <- m-sd
> gppnew <- gpp-v
> plot(gppnew)
> gppnew[gppnew>0] <- 3
> gppnew[gppnew<0] <- 1
> gppnew[gppnew==0] <- 2
> plot(gppnew)
> writeRaster(gppnew, filename=paste("L:/NDVI_WGS84/gppnew_mean_sd",sep=''),
> format="GTiff", overwrite=TRUE)
>
> #Option2_percentile
> max<- cellStats(gpp, 'max')
> gppnew2 <-(gpp/max)*100
> min <- cellStats(gpp, 'min')
> m <- matrix(c(min, 15, 1,
>              15, 85, 2,
>              85, Inf, 3), ncol=3, byrow=TRUE)
> gppnew_percentile <- reclassify(gppnew2, m)
>
> hist(as.matrix(gppnew_percentile))
> writeRaster(gppnew_percentile,
> filename=paste("L:/NDVI_WGS84/gppnew_percentile",sep=''),
> format="GTiff", overwrite=TRUE)
> plot(gppnew_percentile)
>
> ​Thanks for your help
>
> [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Ben Tupper
Bigelow Laboratory for Ocean Sciences
60 Bigelow Drive, P.O. Box 380
East Boothbay, Maine 04544
http://www.bigelow.org

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: raster - unrotate?

Thu, 06/22/2017 - 06:28
Hi,

Wow!  A treasure from the past!

Hmmm.   I see what you mean; it might be best to snip-clip-and-smush from the native presentation.  We are testing out the performance of the dineof function in the sinkr package (https://github.com/marchtaylor/sinkr/blob/master/R/dineof.R) to interpolate patchy chlorophyl data.  

Thanks for the tip!

Cheers,
Ben

> On Jun 22, 2017, at 4:46 AM, Michael Sumner <[hidden email]> wrote:
>
>
> It used to do the inverse, and I preserved a copy of the old behaviour here:
>
> https://github.com/AustralianAntarcticDivision/raadtools/blob/master/R/utils.R#L9 <https://github.com/AustralianAntarcticDivision/raadtools/blob/master/R/utils.R#L9>
>
> You're probably best to learn from how it approaches it rather than just use it, since it's not a general capability, just works for those very specific cases.
>
> I just tested it and it still seems to work fine, I used
>
> oceandata.sci.gsfc.nasa.gov/MODISA/Mapped/Monthly/9km/chlor/A20021822002212.L3m_MO_CHL_chlor_a_9km.nc <http://oceandata.sci.gsfc.nasa.gov/MODISA/Mapped/Monthly/9km/chlor/A20021822002212.L3m_MO_CHL_chlor_a_9km.nc>
>
> obtainable with
>
> https://oceandata.sci.gsfc.nasa.gov/cgi/getfile/A20021822002212.L3m_MO_CHL_chlor_a_9km.nc <https://oceandata.sci.gsfc.nasa.gov/cgi/getfile/A20021822002212.L3m_MO_CHL_chlor_a_9km.nc>
>
> If you're extracting points you might as well just convert them into the "Atlantic" range, but it's painful to get that right for polygons or lines (and very hard to generalize once you get it working anyway).
>
> For simple regions like the one you show I'd probably split it into two extents() and use those - but depends on your workflow, all these options are useful here and there.
>
> Cheers, Mike.
>
>
> On Thu, 22 Jun 2017 at 18:30 Ben Tupper <[hidden email] <mailto:[hidden email]>> wrote:
> Hello,
>
> We have rasters that span [-180, 180] from NASA's Ocean Color (https://oceancolor.gsfc.nasa.gov/ <https://oceancolor.gsfc.nasa.gov/>)  datasets.  We are trying to extract a region that spans 100E to 90W, that is 100 to -90.  The region 'wraps' across the edges as shown by the plot at the address below.
>
> https://dl.dropboxusercontent.com/u/8433654/sst.png <https://dl.dropboxusercontent.com/u/8433654/sst.png>
>
>
> library(raster)
> uri <- 'https://oceandata.sci.gsfc.nasa.gov:443/opendap/MODISA/L3SMI/2016/001/A20160012016032.L3m_R32_SST_sst_9km.nc <https://oceandata.sci.gsfc.nasa.gov/opendap/MODISA/L3SMI/2016/001/A20160012016032.L3m_R32_SST_sst_9km.nc>'
> R <- raster::raster(uri, varname = 'sst')
>
> R
> #class       : RasterLayer
> #dimensions  : 180, 360, 64800  (nrow, ncol, ncell)
> #resolution  : 1, 1  (x, y)
> #extent      : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
> #coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
> #data source : in memory
> #names       : layer
> #values      : 1.482572e-05, 0.9999928  (min, max)
>
> plot(R)
> lines(c(180, 100, 100, 180), c(80,80,0,0))
> lines(c(-180,-90,-90,-180), c(80,80,0,0))
>
> I see that there is the nice raster::rotate() function to rotate raster coordinates from [0,360] to [-180,180]. That would make extracting the region super easy if the inverse were available.  Is there an equivalent way to implement the inverse or raster::rotate()?  I could dig into the source of raster::rotate() to see how to code one up, but I hate like heck to reinvent the wheel.
>
> Thanks!
> Ben
>
>
> Ben Tupper
> Bigelow Laboratory for Ocean Sciences
> 60 Bigelow Drive, P.O. Box 380
> East Boothbay, Maine 04544
> http://www.bigelow.org <http://www.bigelow.org/>
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email] <mailto:[hidden email]>
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo <https://stat.ethz.ch/mailman/listinfo/r-sig-geo>
> --
> Dr. Michael Sumner
> Software and Database Engineer
> Australian Antarctic Division
> 203 Channel Highway
> Kingston Tasmania 7050 Australia
>
Ben Tupper
Bigelow Laboratory for Ocean Sciences
60 Bigelow Drive, P.O. Box 380
East Boothbay, Maine 04544
http://www.bigelow.org




        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: sf 0.5-0 binary for macOS SIerra?

Thu, 06/22/2017 - 06:02
Please also follow:

https://github.com/edzer/sfr/issues/393

and consider trying:

http://www.kyngchaos.com/software/frameworks

for a source sf install until things stabilise.

Roger

On Thu, 22 Jun 2017, Kent Johnson wrote:

> Is there a simple features 0.5-0 binary for macOS Sierra (10.12)? I am
> having trouble installing it from source. CRAN shows
> OS X El Capitan binaries:r-release: sf_0.4-3.tgz
> OS X Mavericks binaries:r-oldrel: sf_0.5-0.tgz
>
> Can I use the Mavericks binary?
>
> Thanks,
> Kent
>
> [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
--
Roger Bivand
Department of Economics, Norwegian School of Economics,
Helleveien 30, N-5045 Bergen, Norway.
voice: +47 55 95 93 55; e-mail: [hidden email]
Editor-in-Chief of The R Journal, https://journal.r-project.org/index.html
http://orcid.org/0000-0003-2392-6140
https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Roger Bivand
Department of Economics
Norwegian School of Economics
Helleveien 30
N-5045 Bergen, Norway

sf 0.5-0 binary for macOS SIerra?

Thu, 06/22/2017 - 05:54
Is there a simple features 0.5-0 binary for macOS Sierra (10.12)? I am
having trouble installing it from source. CRAN shows
OS X El Capitan binaries:r-release: sf_0.4-3.tgz
OS X Mavericks binaries:r-oldrel: sf_0.5-0.tgz

Can I use the Mavericks binary?

Thanks,
Kent

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: raster - unrotate?

Thu, 06/22/2017 - 03:46
It used to do the inverse, and I preserved a copy of the old behaviour
here:

https://github.com/AustralianAntarcticDivision/raadtools/blob/master/R/utils.R#L9

You're probably best to learn from how it approaches it rather than just
use it, since it's not a general capability, just works for those very
specific cases.

I just tested it and it still seems to work fine, I used

oceandata.sci.gsfc.nasa.gov/MODISA/Mapped/Monthly/9km/chlor/A20021822002212.L3m_MO_CHL_chlor_a_9km.nc


obtainable with

https://oceandata.sci.gsfc.nasa.gov/cgi/getfile/A20021822002212.L3m_MO_CHL_chlor_a_9km.nc


If you're extracting points you might as well just convert them into the
"Atlantic" range, but it's painful to get that right for polygons or lines
(and very hard to generalize once you get it working anyway).

For simple regions like the one you show I'd probably split it into two
extents() and use those - but depends on your workflow, all these options
are useful here and there.

Cheers, Mike.


On Thu, 22 Jun 2017 at 18:30 Ben Tupper <[hidden email]> wrote:

> Hello,
>
> We have rasters that span [-180, 180] from NASA's Ocean Color (
> https://oceancolor.gsfc.nasa.gov/)  datasets.  We are trying to extract a
> region that spans 100E to 90W, that is 100 to -90.  The region 'wraps'
> across the edges as shown by the plot at the address below.
>
> https://dl.dropboxusercontent.com/u/8433654/sst.png
>
>
> library(raster)
> uri <- '
> https://oceandata.sci.gsfc.nasa.gov:443/opendap/MODISA/L3SMI/2016/001/A20160012016032.L3m_R32_SST_sst_9km.nc
> '
> R <- raster::raster(uri, varname = 'sst')
>
> R
> #class       : RasterLayer
> #dimensions  : 180, 360, 64800  (nrow, ncol, ncell)
> #resolution  : 1, 1  (x, y)
> #extent      : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
> #coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
> #data source : in memory
> #names       : layer
> #values      : 1.482572e-05, 0.9999928  (min, max)
>
> plot(R)
> lines(c(180, 100, 100, 180), c(80,80,0,0))
> lines(c(-180,-90,-90,-180), c(80,80,0,0))
>
> I see that there is the nice raster::rotate() function to rotate raster
> coordinates from [0,360] to [-180,180]. That would make extracting the
> region super easy if the inverse were available.  Is there an equivalent
> way to implement the inverse or raster::rotate()?  I could dig into the
> source of raster::rotate() to see how to code one up, but I hate like heck
> to reinvent the wheel.
>
> Thanks!
> Ben
>
>
> Ben Tupper
> Bigelow Laboratory for Ocean Sciences
> 60 Bigelow Drive, P.O. Box 380
> East Boothbay, Maine 04544
> http://www.bigelow.org
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> --
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Thresholds & reclassify raster

Thu, 06/22/2017 - 03:04
​Dear all​
,
I am requesting for your help to reclassify the gpp.tif (attached) raster
in R. Not sure my script here is doing the right thing.

​#Data structure​

class       : RasterBrick
dimensions  : 1243, 1409, 1751387, 1  (nrow, ncol, ncell, nlayers)
resolution  : 0.00225804, 0.00225804  (x, y)
extent      : -25.53282, -22.35125, 14.54232, 17.34906  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
data source : L:\NDVI_WGS84\gpp.tif
names       :        gpp
min values  : -0.2176632
max values  :   2.989293


​###​

​I want to ​
use the statistical distribution of gpp.tif values set threshold for change
assignment. The thresholds should be set by the following two
​ ​
options:

(Option1) => use  mean and StD and assign all pixel values  < mean - StD as
negative change (class 1) and all values  > mean + StD as positive
change (class
3); all values in between would be assigned 'neutral' for the final
tabulation (class2).

(option2) => what may come to almost the same result is to do the
cumulative histogram of the calculated GPPproxy pixels and assign all
pixels below the 15 percentile to negative change (class1) and all pixel
above the 85 percentile to positive change (class3) and again the range in
between (15 to 85 percentile) would be assigned neutral (class2).

#R script

# classification
#Option1_mean_sd
gpp <- brick('L:/NDVI_WGS84/gpp.tif')
hist(as.matrix(GPP))
plot(GPPproxy)
m <- cellStats(gpp, 'mean')
sd <- cellStats(gpp, 'sd')
v <- m-sd
gppnew <- gpp-v
plot(gppnew)
gppnew[gppnew>0] <- 3
gppnew[gppnew<0] <- 1
gppnew[gppnew==0] <- 2
plot(gppnew)
writeRaster(gppnew, filename=paste("L:/NDVI_WGS84/gppnew_mean_sd",sep=''),
format="GTiff", overwrite=TRUE)

#Option2_percentile
max<- cellStats(gpp, 'max')
gppnew2 <-(gpp/max)*100
min <- cellStats(gpp, 'min')
m <- matrix(c(min, 15, 1,
              15, 85, 2,
              85, Inf, 3), ncol=3, byrow=TRUE)
gppnew_percentile <- reclassify(gppnew2, m)

hist(as.matrix(gppnew_percentile))
writeRaster(gppnew_percentile,
filename=paste("L:/NDVI_WGS84/gppnew_percentile",sep=''),
format="GTiff", overwrite=TRUE)
plot(gppnew_percentile)

​Thanks for your help

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

raster - unrotate?

Wed, 06/21/2017 - 10:41
Hello,

We have rasters that span [-180, 180] from NASA's Ocean Color (https://oceancolor.gsfc.nasa.gov/)  datasets.  We are trying to extract a region that spans 100E to 90W, that is 100 to -90.  The region 'wraps' across the edges as shown by the plot at the address below.

https://dl.dropboxusercontent.com/u/8433654/sst.png


library(raster)
uri <- 'https://oceandata.sci.gsfc.nasa.gov:443/opendap/MODISA/L3SMI/2016/001/A20160012016032.L3m_R32_SST_sst_9km.nc'
R <- raster::raster(uri, varname = 'sst')

R
#class       : RasterLayer
#dimensions  : 180, 360, 64800  (nrow, ncol, ncell)
#resolution  : 1, 1  (x, y)
#extent      : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
#coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
#data source : in memory
#names       : layer
#values      : 1.482572e-05, 0.9999928  (min, max)

plot(R)
lines(c(180, 100, 100, 180), c(80,80,0,0))
lines(c(-180,-90,-90,-180), c(80,80,0,0))

I see that there is the nice raster::rotate() function to rotate raster coordinates from [0,360] to [-180,180]. That would make extracting the region super easy if the inverse were available.  Is there an equivalent way to implement the inverse or raster::rotate()?  I could dig into the source of raster::rotate() to see how to code one up, but I hate like heck to reinvent the wheel.

Thanks!
Ben


Ben Tupper
Bigelow Laboratory for Ocean Sciences
60 Bigelow Drive, P.O. Box 380
East Boothbay, Maine 04544
http://www.bigelow.org

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: Kriging prediction variances on gstat

Wed, 06/21/2017 - 06:20
Thanks!


> Il giorno 21 giu 2017, alle ore 10:52, Edzer Pebesma <[hidden email]> ha scritto:
>
> citation("gstat")
>
> points you to the literature.
>
>> On 21/06/17 10:16, Carlo Cavalieri wrote:
>> Hi all,
>>
>> I would like to know how the kriging variance is computed by gstat, because on my analysis I need to use a cell-specific covariate depended scale parameter rather than a unique one.
>> If the overall scale of the origins variance is fixed by the sill of the theoretical variogram model used to make the interpolation I can divide the estimated kriging variance by it and then multiply by the cell specific value, but I wanted to make sure that the sill value (partial sill + nugget) is used.
>>
>> Thank you
>> Carlo
>>    [[alternative HTML version deleted]]
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> [hidden email]
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>
> --
> Edzer Pebesma
> Institute for Geoinformatics  (ifgi),  University of Münster
> Heisenbergstraße 2, 48149 Münster, Germany; +49 251 83 33081
> Journal of Statistical Software:   http://www.jstatsoft.org/
> Computers & Geosciences:   http://elsevier.com/locate/cageo/
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: Kriging prediction variances on gstat

Wed, 06/21/2017 - 03:52
citation("gstat")

points you to the literature.

On 21/06/17 10:16, Carlo Cavalieri wrote:
> Hi all,
>
> I would like to know how the kriging variance is computed by gstat, because on my analysis I need to use a cell-specific covariate depended scale parameter rather than a unique one.
> If the overall scale of the origins variance is fixed by the sill of the theoretical variogram model used to make the interpolation I can divide the estimated kriging variance by it and then multiply by the cell specific value, but I wanted to make sure that the sill value (partial sill + nugget) is used.
>
> Thank you
> Carlo
> [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> --
Edzer Pebesma
Institute for Geoinformatics  (ifgi),  University of Münster
Heisenbergstraße 2, 48149 Münster, Germany; +49 251 83 33081
Journal of Statistical Software:   http://www.jstatsoft.org/
Computers & Geosciences:   http://elsevier.com/locate/cageo/


_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
signature.asc (484 bytes) Download Attachment

Re: Issues installing rgeos on linux server

Wed, 06/21/2017 - 03:31
Would it be helpful if we kept a docker file somewhere, that shows how
to install rgeos/rgdal with custom installed libraries? For sf, I
maintain one here:

https://github.com/edzer/sfr/tree/master/inst/docker/custom

On 20/06/17 22:35, Matthew Simpson wrote:
> Thank you, this was very helpful. It turned out the issue was that I was
> using the development version - I had cloned the github repository and
> installed from there. Downloading the tarball for the latest stable release
> solved the problem. I went ahead and did the same thing for GDAL and proj.4
> and reinstalled rgdal just to be safe. And sure enough, before reinstalling
> rgdal was complaining about missing proj_defs.dat, even though it existed
> in the proj.4/share folder. But now it's working just fine.
>
> For posterity, my full commands to get GEOS, GDAL, proj.4, rgdal, and rgeos
> installed from source on a linux server (CentOS 7.3.1611) are below.
>
> cd /path/to/src/geos/
> wget http://download.osgeo.org/geos/geos-3.6.1.tar.bz2
> tar -vxjf geos-3.6.1.tar.bz2
> cd geos-3.6.1
> ./configure --prefix=/path/to/install/geos
> make
> make check
> make install
>
> cd /path/to/src/gdal
> wget http://download.osgeo.org/gdal/2.2.0/gdal-2.2.0.tar.gz
> tar -xvzf gdal-2.2.0.tar.gz
> cd gdal-2.2.0
> ./autogen.sh
> ./configure --prefix=/path/to/install/gdal
> make
> ## there is no make check
> make install
>
> cd /path/to/src/proj.4/
> wget http://download.osgeo.org/proj/proj-4.9.3.tar.gz
> tar -xvzf proj-4.9.3.tar.gz
> cd proj-4.9.3
> ./configure --prefix=/path/to/install/proj.4
> make
> make check
> make install
>
> echo 'export
> PATH=/path/to/install/gdal/bin:/path/to/install/proj.4/bin:/path/to/install/geos/bin:$PATH'
>>> ~/.bash_profile
> echo 'export
> LD_LIBRARY_PATH=/path/to/install/gdal/lib:/path/to/install/proj.4/lib:/path/to/install/geos/lib:$LD_LIBRARY_PATH'
>>> ~/.bash_profile
> . ~/.bash_profile
>
> ## in R run:
> ## need this installed from source to prevent some compiler issues that
> would cause R to segfault if rgdal was loaded and you tried to install
> certain packages
> install.packages("Rcpp", type = "source", repos = "
> https://cloud.r-project.org")
>
> install.packages("rgdal", type = "source", repos = "
> https://cloud.r-project.org",
> configure.args = "--with-proj-include=/path/to/install/proj.4/include
> --with-proj-lib=/path/to/install/proj.4/lib
> --with-proj-share=/path/to/install/proj.4/share/proj
> --with-gdal-config=/path/to/install/gdal/bin/gdal-config")
>
> install.packages("rgeos", type = "source", repos = "
> https://cloud.r-project.org",
>          configure.args =
> "--with-geos-config=/path/to/install/geos/bin/geos-config")
>
>
> Note that --with-proj-share=/path/to/install/proj.4/share did NOT work.
>
>
> On Tue, Jun 20, 2017 at 4:26 AM, Roger Bivand <[hidden email]> wrote:
>
>> On Tue, 20 Jun 2017, Matthew Simpson wrote:
>>
>> I'm installing some spatial packages on a linux server and due to
>>> university constraints, cannot install anything using a package manager.
>>> I've managed to get rgdal installed successfully, but rgeos is causing
>>> problems.
>>>
>>> I've installed the geos library already, but when I attempt to install
>>> rgeos within R, it complains that that it cannot run C compiled programs.
>>>
>>
>> Could you please say how you installed GEOS? Could you try out geos-config
>> with each of the settings? Do you need 3.7.0dev (I haven't tried this as a
>> source of difficulty as it seems implausible, but neither rgeos nor sf need
>> bleeding edge functionality from GEOS)? Could you try to install sf from
>> source (it also uses GEOS)? Is this the GEOS geos-config -m library bug
>> (probably not)? Unpacking rgeos, running ./configure in the package root, I
>> see (Fedora 25, R 3.4.0):
>>
>> $ ./configure
>> configure: CC: gcc
>> configure: CXX: g++
>> configure: rgeos: 0.3-23
>> checking for /usr/bin/svnversion... yes
>> configure: svn revision: 547
>> checking for geos-config... /usr/local/bin/geos-config
>> checking geos-config usability... yes
>> configure: GEOS version: 3.6.1
>> checking geos version at least 3.2.0... yes
>> checking geos-config clibs... yes
>> checking for gcc... gcc
>> checking whether the C compiler works... yes
>> checking for C compiler default output file name... a.out
>> checking for suffix of executables...
>> checking whether we are cross compiling... no
>> checking for suffix of object files... o
>> checking whether we are using the GNU C compiler... yes
>> checking whether gcc accepts -g... yes
>> checking for gcc option to accept ISO C89... none needed
>> checking how to run the C preprocessor... gcc -E
>> checking for grep that handles long lines and -e... /usr/bin/grep
>> checking for egrep... /usr/bin/grep -E
>> checking for ANSI C header files... yes
>> checking for sys/types.h... yes
>> checking for sys/stat.h... yes
>> checking for stdlib.h... yes
>> checking for string.h... yes
>> checking for memory.h... yes
>> checking for strings.h... yes
>> checking for inttypes.h... yes
>> checking for stdint.h... yes
>> checking for unistd.h... yes
>> checking geos_c.h usability... yes
>> checking geos_c.h presence... yes
>> checking for geos_c.h... yes
>> checking geos: linking with libgeos_c... yes
>> configure: PKG_CPPFLAGS:  -I/usr/local/include
>> configure: PKG_LIBS:  -L/usr/local/lib -lgeos -L/usr/local/lib -lgeos_c
>> configure: creating ./config.status
>> config.status: creating src/Makevars
>>
>> but you fail about half-way down at the standard croos-compilation test.
>> For rgdal/configure - you said your rgdal install on the same system
>> succeeded, I see:
>>
>> ...
>> checking for gcc... gcc
>> checking whether the C compiler works... yes
>> checking for C compiler default output file name... a.out
>> checking for suffix of executables...
>> checking whether we are cross compiling... no
>> checking for suffix of object files... o
>> checking whether we are using the GNU C compiler... yes
>> checking whether gcc accepts -g... yes
>> checking for gcc option to accept ISO C89... none needed
>> checking how to run the C preprocessor... gcc -E
>> ...
>>
>> What do you see (just running ./configure in the package root)?
>>
>> Roger
>>
>>
>>> I've also tried downloading the tarball myself and running R CMD check,
>>> and
>>> it gives me the same complaint.
>>>
>>> I've successfully reinstalled Rcpp from source, but that seemed to do
>>> nothing.
>>>
>>> Below is: 1) the console input from the attempted installation, 2) shell
>>> output from running R CMD check, and 3) my sessionInfo()
>>>
>>> Any help would be appreciated.
>>>
>>> Matt
>>>
>>> 1) Here is the result of attempted installation: (note: removing type =
>>> "source" changes nothing):
>>>
>>> install.packages("rgeos", type = "source",
>>>>
>>> +         configure.args =
>>> "--with-geos-config=/group/stsn/shared_libs/geos/bin/geos-config",
>>> +         repos = "https://cloud.r-project.org")
>>> Installing package into ‘/home/msimpson/R/x86_64-pc-li
>>> nux-gnu-library/3.3’
>>> (as ‘lib’ is unspecified)
>>> trying URL 'https://cloud.r-project.org/src/contrib/rgeos_0.3-23.tar.gz'
>>> Content type 'application/x-gzip' length 257486 bytes (251 KB)
>>> ==================================================
>>> downloaded 251 KB
>>>
>>> * installing *source* package ‘rgeos’ ...
>>> ** package ‘rgeos’ successfully unpacked and MD5 sums checked
>>> configure: CC: /usr/bin/gcc -std=gnu99
>>> configure: CXX: /usr/bin/g++
>>> configure: rgeos: 0.3-23
>>> checking for /usr/bin/svnversion... no
>>> configure: svn revision: 546
>>> configure: geos-config set to /group/stsn/shared_libs/geos/b
>>> in/geos-config
>>> checking geos-config exists... yes
>>> checking geos-config executable... yes
>>> checking geos-config usability... yes
>>> configure: GEOS version: 3.7.0dev
>>> checking geos version at least 3.2.0... yes
>>> checking geos-config clibs... yes
>>> checking for gcc... /usr/bin/gcc -std=gnu99
>>> checking whether the C compiler works... yes
>>> checking for C compiler default output file name... a.out
>>> checking for suffix of executables...
>>> checking whether we are cross compiling... configure: error: in
>>> `/tmp/RtmppP5J2b/R.INSTALL236bb22fea27c/rgeos':
>>> configure: error: cannot run C compiled programs.
>>> If you meant to cross compile, use `--host'.
>>> See `config.log' for more details
>>> ERROR: configuration failed for package ‘rgeos’
>>> * removing ‘/home/msimpson/R/x86_64-pc-linux-gnu-library/3.3/rgeos’
>>>
>>> The downloaded source packages are in
>>>        ‘/tmp/RtmpwkU04N/downloaded_packages’
>>> Warning message:
>>> In install.packages("rgeos", type = "source", configure.args =
>>> "--with-geos-config=/group/stsn/shared_libs/geos/bin/geos-config",  :
>>>  installation of package ‘rgeos’ had non-zero exit status
>>>
>>> 2) Here's the result from running R CMD check:
>>>
>>> R CMD check
>>> --install-args='--configure-args=--with-geos-config=/group/
>>> stsn/shared_libs/geos/bin/geos-config'
>>> rgeos_0.3-23.tar.gz
>>> * using log directory ‘/group/stsn/src/test/rgeos.Rcheck’
>>> * using R version 3.3.1 (2016-06-21)
>>> * using platform: x86_64-pc-linux-gnu (64-bit)
>>> * using session charset: UTF-8
>>> * checking for file ‘rgeos/DESCRIPTION’ ... OK
>>> * this is package ‘rgeos’ version ‘0.3-23’
>>> * checking package namespace information ... OK
>>> * checking package dependencies ... OK
>>> * checking if this is a source package ... OK
>>> * checking if there is a namespace ... OK
>>> * checking for executable files ... OK
>>> * checking for hidden files and directories ... OK
>>> * checking for portable file names ... OK
>>> * checking for sufficient/correct file permissions ... OK
>>> * checking whether package ‘rgeos’ can be installed ... ERROR
>>> Installation failed.
>>> See ‘/group/stsn/src/test/rgeos.Rcheck/00install.out’ for details.
>>> * DONE
>>>
>>> Status: 1 ERROR
>>> See
>>>  ‘/group/stsn/src/test/rgeos.Rcheck/00check.log’
>>> for details.
>>>
>>> cat /group/stsn/src/test/rgeos.Rcheck/00install.out
>>> * installing *source* package ‘rgeos’ ...
>>> ** package ‘rgeos’ successfully unpacked and MD5 sums checked
>>> configure: CC: /usr/bin/gcc -std=gnu99
>>> configure: CXX: /usr/bin/g++
>>> configure: rgeos: 0.3-23
>>> checking for /usr/bin/svnversion... no
>>> configure: svn revision: 546
>>> configure: geos-config set to /group/stsn/shared_libs/geos/b
>>> in/geos-config
>>> checking geos-config exists... yes
>>> checking geos-config executable... yes
>>> checking geos-config usability... yes
>>> configure: GEOS version: 3.7.0dev
>>> checking geos version at least 3.2.0... yes
>>> checking geos-config clibs... yes
>>> checking for gcc... /usr/bin/gcc -std=gnu99
>>> checking whether the C compiler works... yes
>>> checking for C compiler default output file name... a.out
>>> checking for suffix of executables...
>>> checking whether we are cross compiling... configure: error: in
>>> `/group/stsn/src/test/rgeos.Rcheck/00_pkg_src/rgeos':
>>> configure: error: cannot run C compiled programs.
>>> If you meant to cross compile, use `--host'.
>>> See `config.log' for more details
>>> ERROR: configuration failed for package ‘rgeos’
>>> * removing ‘/group/stsn/src/test/rgeos.Rcheck/rgeos’
>>>
>>> 3) And here's my session info:
>>>
>>> sessionInfo()
>>> R version 3.3.1 (2016-06-21)
>>> Platform: x86_64-pc-linux-gnu (64-bit)
>>> Running under: CentOS Linux 7 (Core)
>>>
>>> locale:
>>> [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C
>>> [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8
>>> [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8
>>> [7] LC_PAPER=en_US.UTF-8       LC_NAME=C
>>> [9] LC_ADDRESS=C               LC_TELEPHONE=C
>>> [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
>>>
>>> attached base packages:
>>> [1] stats     graphics  grDevices utils     datasets  methods   base
>>>
>>>         [[alternative HTML version deleted]]
>>>
>>> _______________________________________________
>>> R-sig-Geo mailing list
>>> [hidden email]
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>
>>
>> --
>> Roger Bivand
>> Department of Economics, Norwegian School of Economics,
>> Helleveien 30, N-5045 Bergen, Norway.
>> voice: +47 55 95 93 55; e-mail: [hidden email]
>> Editor-in-Chief of The R Journal, https://journal.r-project.org/index.html
>> http://orcid.org/0000-0003-2392-6140
>> https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en
>
> [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> --
Edzer Pebesma
Institute for Geoinformatics  (ifgi),  University of Münster
Heisenbergstraße 2, 48149 Münster, Germany; +49 251 83 33081
Journal of Statistical Software:   http://www.jstatsoft.org/
Computers & Geosciences:   http://elsevier.com/locate/cageo/


_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
signature.asc (484 bytes) Download Attachment

Kriging prediction variances on gstat

Wed, 06/21/2017 - 03:16
Hi all,

I would like to know how the kriging variance is computed by gstat, because on my analysis I need to use a cell-specific covariate depended scale parameter rather than a unique one.
If the overall scale of the origins variance is fixed by the sill of the theoretical variogram model used to make the interpolation I can divide the estimated kriging variance by it and then multiply by the cell specific value, but I wanted to make sure that the sill value (partial sill + nugget) is used.

Thank you
Carlo
        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Problems with using dismo package on HPC

Tue, 06/20/2017 - 16:30
Hi all,

I am trying to use the dismo package on an HPC and am running into some problems. The code that I used worked fine until a few days ago, then I started getting the error message

Error in rgdal::getRasterData(con, offset = offs, region.dim = c(1, nc),  :
Failure during raster IO

Not sure how to proceed on this. The code works on my desktop computer, and like I said, worked until a few days ago on the HPC. I have read/write privileges to the outfolder location,  all the environmental layers are the same extent and resolution, and the spatial points are in the same coordinate system.

Any help would be greatly appreciated.
Wade

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Pages