Tag Info

New answers tagged

1

I kind misunderstand your question. but anyway... First of all you need to download the ggmap library in order to work with it later. Try this code out ! And you can use Google API in order to plot GPS coordinate in a map!


0

Here I post my answer: I missed out how to convert SpatialPolygon to SpatialPolygonDataFrame and difference between them. SpatialPolygon doesn't contain attribute table, SpatialPolygonDataFrame does. So attribute table can be easily add from my points attribute table. No need to use over function :-D # re-written code library(sp) library(deldir) ...


4

Here are some ideas. With base plot you can do plot(x, interpolate=TRUE) You can also resample your data y <- disaggregate(x, 5, method='bilinear') Or indeed smooth it using a focal operation y <- focal(x, w=matrix(1, 5, 5), mean) Or a combination y <- disaggregate(x, 5) y <- focal(y, w=matrix(1, 5, 5), mean) The question ...


2

You cannot have R objects called "2000", so presumably these are fake names? Your example actually should work, so you may want to double check why you think that the results are incorrect. @aaryno's approach should work. I would do this: library(raster) s <- stack(r2000, r2001, r2002, r2003, r2004, r2005) x <- reclassify(s, cbind(0, NA)) r <- ...


2

I wanted to point out that you can rewrite the mean function, mean, which you can write yourself to do anything you want, including ditch the 0 values and calculate the mean. For example, if you want to ignore 0s: meanIgnoringZeroes <- function(x) { mean(x[x!=0],na.rm=T) } Then you can pass the function, meanIgnoringZeroes to overlay: mean <- ...


-2

Perhaps you should try this : new_file=na.omit(file). I think with this you can ignore NAs.


1

Although this may not be directly apparent from the documentation, the third (unnamed) argument to your krige.cv call is assumed to be a variogram model, and is absorbed by the argument model. Since dec_vor is the outcome of a call to krige, it is not, hence the error message. It is not clear to me why you would pass dec_vor to krige.cv, as krige.cv loops ...


0

Converts my solution to a comment. Maybe more words makes it less trivial? Regardless, all these answers may work, but the SuperOverlay format is horrible, and the quads thing is pretty limiting/crude. I reverse engineered an output from OKMap... And you could use that, but I posted a script for ArcGIS here:Exporting 3GB ArcGIS Raster to KML without losing ...


4

Compute focal means of the indicators of each vegetation type. At each cell, these give the proportions of the types. Multiply each by its negative logarithm and sum: that's the diversity index. You will find that even for large numbers of categories (even into the hundreds), this is fantastically faster than the brute-force method of tabulating each ...


5

Here's a simple implementation for you. Edit: fixed focal weights matrix to exclude 0s as per whuber's comments in his answer. library(raster) # Example Data set.seed(1) r <- raster(matrix(sample(1:10, 100, replace=T), 10, 10)) # Calculate a weights matrix, and reset elements to 0s and 1s # rather than true weights fw <- focalWeight(r, 0.2, ...


1

OK, here comes the correct answer: Make sure that rgdal (version >= 1.0.4) is installed install.packages('rgdal') packageVersion('rgdal') [1] ‘1.0.4’ Make sure that gdal (version >= 1.11.0) is installed library(rgdal) getGDALVersionInfo() [1] "GDAL 1.11.2, released 2015/02/10" Make sure that gdal is compiled with Expat/OSM and SQLite support: ...


1

You need to convert your SpatialPolygons class to a SpatialPolygonsDataFrame class. For example: require(rgdal) require(rgeos) # Read shapefile shp = 'C:/temp/myshp.shp' myshp = readOGR(shp, layer = basename(strsplit(shp, "\\.")[[1]])[1]) # Read shapefile attributes df = data.frame(myshp) # Simplify geometry using rgeos simplified = gSimplify(myshp, tol ...


1

Coerce your object to the appropriate Spatial*DataFrame-class (Points/Lines/Polygons), e.g. for SpatialPolygons using as(x, "SpatialPolygonsDataFrame" ): R> l <- readWKT("LINESTRING(0 7,1 6,2 1,3 4,4 1,5 7,6 6,7 4,8 6,9 4)") R> x1 <- gSimplify(p, tol=10) R> class(x1) [1] "SpatialPolygons" attr(,"package") [1] "sp" R> x2 <- as(x, ...


0

I got the solution. For creating a new grid we first need the bounding box and the pixel size to which it is to be changed. The below code would form a new grid of any resolution: #defining spatial grid data frame S=5 #scale to which it is to be converted psize<-psdeg[1]/S # psdeg is the pixel size of old grid bb<-c@bbox #bounding box where "c" ...


0

You can do: library(raster) d <- raster(c) dd <- disaggregate(d, 5)


0

You can use rgeos and raster functions to accomplish this #### set up example data library(raster) library(rgeos) p1 <- rbind(c(-180,-20), c(-140,55), c(10, 0), c(-140,-60), c(-180,-20)) hole <- rbind(c(-150,-20), c(-100,-10), c(-110,20), c(-150,-20)) p1 <- list(p1, hole) p2 <- rbind(c(-10,0), c(140,60), c(160,0), c(140,-55), c(-10,0)) p3 <- ...


0

Your data are invalid: > bbox(worldmapLines) min max x -179.95721 190.29080 y -85.44308 83.57391 and the plot is correct, given that you provided it with wrong data. With help from Roger Bivand: a corrected version of this dataset is available in maptools, try: library(maptools) data(wrld_simpl) plot(spTransform(wrld_simpl, ...


2

Your code has too many brackets for one thing. I count three ( and two ). Fix that, and using file.choose() should work if the user chooses the .shp file. Note that: file.choose() only produces a dialog on Windows (I think) - on Linux its a text prompt. I don't know if it works in RStudio... The maptools functions readShape* usually fail to read ...


1

Got it working using the fmi package in R. https://github.com/rOpenGov/fmi/blob/master/vignettes/fmi_tutorial.md


2

You can use rgdal to access feature classes in Esri file geodatabases. require(rgdal) # The input file geodatabase fgdb = "C:/path/to/your/filegeodatabase.gdb" # List all feature classes in a file geodatabase subset(ogrDrivers(), grepl("GDB", name)) fc_list = ogrListLayers(fgdb) print(fc_list) # Read the feature class fc = ...


3

So, based on AndreJ's answer above, I whipped this up. cleanSpLinesForProjection <- function(w){ #We'll have to go through this one, line by line... slot(w, "lines") <- lapply(slot(worldmapLines, "lines"), function(x) { coords <- slot(x, "Lines")[[1]]@coords #get the rows with too large longitude values rIDX <- which(coords ...


3

Your data exceeds the 180°E limit in eastern russia and some pacific island. If you want a good looking map in degrees, you have to cut your source data at 179.9°E/W. See my answer here (though it deals with a pacific centered view) QGIS display world country shape files centered on pacific ocean using Robinson, Miller Cylindrical or other projection


2

To select the files you could parse the numbers, but if the file names are as regular as you say (January is '01') I think you can also use their character representation ( '9' > '8' == TRUE ) library(raster) selectedFilesFun <- function(files, start, end, fun) { b <- basename(files) b <- substr(b, 2, 9) i <- b > start & b < ...


1

Have you tried the fmi R package from rOpenGov? This is an example partly from their tutorial and ran for me just fine: install.packages(c("devtools", "sp", "rgdal", "raster")) # if you haven't library(devtools) install_github("rOpenGov/rwfs") install_github("rOpenGov/fmi") library(fmi) library(sp) apiKey <- "ENTER YOUR API KEY HERE" request <- ...


0

Here is what Mikkel suggested (use of max) library(raster) cell100 <- raster(nr=3, nc=3, vals=c(100,0,0,0,100,0,0,0,100)) cell101 <- raster(nr=3, nc=3, vals=c(0,0,101,0,101,0,101,0,0)) r <- max(cell100, cell101) as.matrix(r) # [,1] [,2] [,3] #[1,] 100 0 101 #[2,] 0 101 0 #[3,] 101 0 100 Another (more complex) approach could ...


2

You can use lapply to return the mean(s) for a specific column present in all of your sp @data slots. Create a list object with 2 sp point objects library(sp) data(meuse) coordinates(meuse) <- ~x+y sp.list <- list() sp.list[[1]] <- meuse[1:78,] sp.list[[2]] <- meuse[79:155,] Now we can return the mean(s) for the cadmium column in the ...


1

You could use a list to save all the krige.cv data.frames. I changed your example to the meuse data to make it reproducible. I'm assuming your data is similarly formatted. require(gstat) require(sp) data(meuse) data(meuse.grid) coordinates(meuse) = ~x+y coordinates(meuse.grid) = ~x+y gridded(meuse.grid) = TRUE power = seq(from = 1, to = 1.6, by = 0.1) ...


4

You can make a paletted raster by assigning a colortable in the legend. If you have a raster called r and a data frame like yours above called ctab, with value and red/green/blue colour values, you can do something like this: > ctable = rep(NA,max(ctab$value)+1) > ctable[ctab$value+1] = rgb(ctab$red,ctab$green,ctab$blue,maxColorValue=255) > ...


1

A simple "ifelse" statement should suffice in evaluating a condition. Here we create two random vectors of [1,2] and apply an ifelse to evaluate the condition of if x = 2 and y = 1 THEN change (1) ELSE no change (0). ( x <- round(runif(10, 1, 2)) ) ( y <- round(runif(10, 1, 2)) ) ifelse( x == 2 & y == 1, 1, 0) Since this is just an ...


5

Index each day with its year and julian day. September 1 is 244, March 31 is 90. Not knowing the format of your files, you can probably figure out a date format to use to parse and turn it into a POSIX date (from http://stackoverflow.com/questions/21414847/convert-a-date-vector-into-julian-day-in-r): tmp <- as.POSIXlt("16Jun10", format = "%d%b%y") ...


0

Convert your map to a SpatialPolygons object and then use gIntersection from the rgeos package. This clips your blue polygons exactly to the state boundaries. Inspired by this answer. library(rgeos) library(maps) library(maptools) mmap <- map('state', regions=c('maryland', 'virginia', 'delaware'), fill=TRUE) IDs <- sapply(strsplit(mmap$names, ":"), ...


1

Thank you very much for your responses! Eventually, is very very easy to work with GRASS GIS through R. After you have created in GRASS the location and mapset in which you wish to work, you can type in the GRASS shell: "rstudio &" "&" Helps for working simultaneously in both GRASS GIS and R. Otherwise, the GRASS shell would switch to R. And that is ...


2

have a look at the examples in sp: library(sp) example(aggregate)


2

Do it via raster, as long as your points really are regular and the values are discrete (for dissolve) this should work: library(raster) ## x is the SpatialPointsDataFrame r <- rasterFromXYZ(x) ## create polygons on unique values p <- rasterToPolygons(r, dissolve = TRUE) Here's an example library(raster) x <- as(raster(volcano%/% 10) * 10, ...


0

My solution was to convert the polygons to raster (snapped to the land cover map and with the same resolution) using a unique polygon ID field for raster value. Then I converted this raster to points, and ended up with nearly 35,000 points, with attributes POINTID and GRID_CODE. POINTID is a unique cell identifier and GRID_CODE matches the original polygon ...


2

SOLVED: To retrieve R scripts "examples", they have to be installed using the "Get R scripts from on-line scripts collection" tool under the R scripts-> Tools item in the Processing Toolbox. From there, select and download, one by one, many example.


0

This is not a full solution, but it is too big for a comment. The main problem that you are facing is that the mosaic function appears to operates on a per-band basis. As such, you can't use the last band to determine how it should mosaic the first band. A dirty work-around would be to by-pass the mosaic function entirely: Determine the combined total ...


0

Jeffrey: > update.packages(repos="cran.us.r-project.org", ask=FALSE) Warning: unable to access index for repository cran.us.r-project.org/src/contrib


5

Here is a quick worked example of setting the GRASS environment, reading an on-disk raster, calculating a focal mean (using r.neighbors) and reading the results back into R. Hopefully this will get you started. if (!require(rgrass7)) stop("rgrass7 PACKAGE MISSING") setwd("D:/TMP") # Working directory # Set on-disk raster variable rname <- ...


0

You would have to take a good hard look at the script to make sure the the code structure is appropriate for atmospheric MODIS data and that it is as simple as calling a different exe. Since this is a different data format, I highly doubt that it is as simple as calling a different processing program, but one can hope. As far as calling the ...


1

There is no simple implementation of a Kernal Density Estimate using weights in R. Most of the advice for KDE's are limited to spatial locations only. You can write a function to project results from the ks package to a grid, but this is not entirely straight forward. My best advice is to leverage existing implementations from a GIS. The best option I have ...


0

Take a look at density() function in spatstat package. Official site has a number of manuals and articles about this package (see Documentation). I would recommend to start with Analysing Spatial Point Patterns in R.


0

As pointed out by @cengel, the ppp function in spatstat package will take care of the work. Please check the following link for detailed information: http://www.inside-r.org/packages/cran/spatstat/docs/as.ppp


0

In ArcGIS, this is a simple field calculation and export if I understand the question. Add the points layer to the map and open the attribute table. Click on the top left button and choose Add Field. Give it a name and choose the appropriate data type (int/float/double). Then right-click the new field heading and choose Field Calculator. In the resulting ...


2

See How to open a Shapefile in R? for how to open your shapefile. Once you have it open, you have a SpatialPointsDataFrame with a data.frame inside that contains all your attributes. You can then perform any operation row-wise on the data.frame. require(rgdal) # for shapefiles, first argument of the read/write/info functions is the # directory location, ...


3

I'm not sure whether your LiDAR data are raw (i.e. a point cloud) or have already been interpolated to a raster DEM, but if it is the former, you'll need to interpolate the data to a raster DEM first. The next thing you need to do is to perform a 'Depression Filling' operation. Most GIS will have a tool for this task. Once you've created your filled DEM ...


1

This is even simpler: x = rbind(x1, x2, x3, makeUniqueIDs = TRUE)


1

I ran into the same problem. It is inconvenient to use ArcMap, in my opinion. Instead, I invoked python from R using the following code. It requires you to have Arcpy, unfortunately. Python script: import os,sys import arcpy from arcpy import env from sys import argv ### This is needed to import variables script, featureClass, inFeatures, outLocation, ...


1

The ELKI version of DBSCAN has full support for geodetic distances. Just set the distance function to LatLngDistanceFunction or LngLatDistanceFunction (depending on your data format), and specify your epsilon radius in meters. ELKI also has R*-tree index acceleration, making this type of clustering very fast. Benchmark it against R, and you will see R lose ...


2

I can't comment or upvote.. So here is my little contribution; I would also go on a combination of rgdal and RPostgreSQL. So, same code as @Guillaume, except with a tryCatch that handle more lines, a pseudo-random table name and the use of an unlogged table for better performance. (Note to myself: we can't use TEMP table, because it's not visible from ...



Top 50 recent answers are included