[R-sig-Geo] writing shapefiles using write.pointShape

modern82376 at mypacks.net modern82376 at mypacks.net
Thu Mar 13 22:10:27 CET 2008


I've been running gwr on some survey data and trying to write output back out to point shapefiles. At first, I thought that the shapefile wasn't being written, but after calling the file object I found it in a local settings/temp folder with a random shapefile name. I thought I had specified all that, but I'm obviously not familiar enough with the syntax to get it right. Would someone mind assisting me on this?  Code follows..

> library(spgwr)
> surveypts <- readShapePoints("C:/temp/NE_survey.shp", proj4string = CRS("+proj=aea +lat_0=96w +lat_1=29.5n +lat_2=45.5n"))
> survey.bw <- gwr.sel(GHDENSITY ~ pca1 + PCA2 + PCA3 + PCA4 + PCA5, data=surveypts)
> fitpts <- readShapePoints("C:/temp/fitpts.shp", proj4string = CRS("+proj=aea +lat_0=96w +lat_1=29.5n +lat_2=45.5n"))
> fitpix <- SpatialPixels(fitpts, proj4string = CRS("+proj=aea +lat_0=96w +lat_1=29.5n +lat_2=45.5n"))
> surv.local <- ggwr(GHDENSITY ~ pca1 + PCA2 + PCA3 + PCA4 + PCA5, data = surveypts, bandwidth = survey.bw, gweight = gwr.Gauss, adapt = NULL, fit.points = fitpix, longlat = FALSE)
> surv.local

Call:
ggwr(formula = GHDENSITY ~ pca1 + PCA2 + PCA3 + PCA4 + PCA5, 
    data = surveypts, bandwidth = survey.bw, gweight = gwr.Gauss, 
    adapt = NULL, fit.points = fitpix, longlat = FALSE)
Kernel function: gwr.Gauss 
Fixed bandwidth: 22738.78 
Summary of GWR coefficient estimates:
                   Min.    1st Qu.     Median    3rd Qu.       Max.   Global
X.Intercept. -4.176e+03 -9.067e+01 -3.565e+00  8.029e+01  4.402e+03 -73.5570
pca1         -2.567e+00 -2.289e-02 -1.725e-03  1.674e-02  3.135e-01  -0.0077
PCA2         -7.874e+01 -6.543e-01 -1.934e-01  2.298e-01  2.016e+01  -0.2667
PCA3         -8.844e+01 -1.236e+00 -1.109e-01  6.008e-01  1.469e+02   0.0202
PCA4         -1.272e+01 -3.509e-01  3.005e-01  1.178e+00  3.386e+02   0.8161
PCA5         -1.368e+02 -1.134e+00 -8.233e-03  9.222e-01  4.384e+01   0.6937

> names(surv.local$SDF)

[1] "sum.w"           "X.Intercept."    "pca1"            "PCA2"            "PCA3"           
[6] "PCA4"            "PCA5"            "dispersion"      "response_resids"

> coords <- coordinates(fitpts)
> pca1 <- cbind(coords, surv.local$SDF$pca1)
> summary(pca1)
   coords.x1         coords.x2            V3           
 Min.   :-670120   Min.   :287533   Min.   :-2.566983  
 1st Qu.:-544120   1st Qu.:375533   1st Qu.:-0.022891  
 Median :-417620   Median :464533   Median :-0.001725  
 Mean   :-417620   Mean   :464533   Mean   :-0.047442  
 3rd Qu.:-291120   3rd Qu.:553533   3rd Qu.: 0.016741  
 Max.   :-165120   Max.   :641533   Max.   : 0.313544 

>class(pca1)
[1] "matrix"
> dim(pca1)
[1] 179630      3
> all(dim(dpca1 <- as.data.frame(pca1)) == c(179630,3))
[1] TRUE
> summary(dpca1)
   coords.x1         coords.x2            V3           
 Min.   :-670120   Min.   :287533   Min.   :-2.566983  
 1st Qu.:-544120   1st Qu.:375533   1st Qu.:-0.022891  
 Median :-417620   Median :464533   Median :-0.001725  
 Mean   :-417620   Mean   :464533   Mean   :-0.047442  
 3rd Qu.:-291120   3rd Qu.:553533   3rd Qu.: 0.016741  
 Max.   :-165120   Max.   :641533   Max.   : 0.313544  
> class(dpca1)
[1] "data.frame"

> B1pca1 <- tempfile("C:\temp")
> write.pointShape(coords, df=dpca1, B1pca1, factor2char=TRUE)
(Where's the shapefile??)

>B1pca1
[1] "C:\\DOCUME~1\\gcook\\LOCALS~1\\Temp\\Rtmpijon4d\\C:\temp1eb26e9"


After calling B1pca1, I find a randomly named shapefile (with .dbf and .shx)extensions in the .../Local Settings/Temp/ folder. I'm trying to name the file "B1pca1.shp" and I want it to go in the C:/Temp folder (for now). I've tried a few permutations of the last two commands, but I'm just not getting it right.

Thanks so much :)

Gericke




More information about the R-sig-Geo mailing list