Diff to HTML by rtfpessoa

Files changed (1) show
{orig → migrated}/sppa_mod.R RENAMED
@@ -16,37 +16,42 @@
16
16
  ###################################################
17
17
  ### code chunk number 8: sppa.Rnw:258-259
18
18
  ###################################################
19
- library(maptools)
20
19
 
21
20
 
22
21
  ###################################################
23
22
  ### code chunk number 9: sppa.Rnw:261-263
24
23
  ###################################################
25
- spjpines <- as(japanesepines, "SpatialPoints")
24
+ library(sf)
25
+ oo <- st_as_sf(rescale(japanesepines))
26
+ library(sp)
27
+ spjpines <- as(as(oo[oo$label=="point",], "Spatial"), "SpatialPoints")
26
28
  summary(spjpines)
29
+ oo <- st_as_sf(japanesepines)
30
+ spjpines1 <- as(as(oo[oo$label=="point",], "Spatial"), "SpatialPoints")
27
31
 
28
32
 
29
33
  ###################################################
30
34
  ### code chunk number 10: sppa.Rnw:275-277
31
35
  ###################################################
32
- spjpines1 <- elide(spjpines, scale=TRUE, unitsq=TRUE)
33
- summary(spjpines1)
34
36
 
35
37
 
36
38
  ###################################################
37
39
  ### code chunk number 11: sppa.Rnw:294-296
38
40
  ###################################################
39
- pppjap <- as(spjpines1, "ppp")
40
- summary(pppjap)
41
41
 
42
42
 
43
43
  ###################################################
44
44
  ### code chunk number 12: sppa.Rnw:300-311
45
45
  ###################################################
46
46
  data(redwoodfull)
47
- spred <- as(redwoodfull, "SpatialPoints")
47
+ oo <- st_as_sf(redwoodfull)
48
+ spred <- as(as(oo[oo$label=="point",], "Spatial"), "SpatialPoints")
49
+ summary(spred)
48
50
  data(cells)
49
- spcells <- as(cells, "SpatialPoints")
51
+ summary(cells)
52
+ oo <- st_as_sf(cells)
53
+ spcells <- as(as(oo[oo$label=="point",], "Spatial"), "SpatialPoints")
54
+ summary(spcells)
50
55
  dpp<-data.frame(rbind(coordinates(spjpines1), coordinates(spred),
51
56
  coordinates(spcells)))
52
57
  njap<-nrow(coordinates(spjpines1))
@@ -66,11 +71,10 @@
66
71
  ###################################################
67
72
  ### code chunk number 15: sppa.Rnw:401-406
68
73
  ###################################################
69
- library(rgdal)
70
- spasthma <- readOGR("spasthma.shp", "spasthma", stringsAsFactors=TRUE)
71
- spbdry <- readOGR("spbdry.shp", "spbdry")
72
- spsrc <- readOGR("spsrc.shp", "spsrc")
73
- sproads <- readOGR("sproads.shp", "sproads")
74
+ spasthma <- as(st_read("spasthma.shp", stringsAsFactors=TRUE), "Spatial")
75
+ spbdry <- as(st_read("spbdry.shp"), "Spatial")
76
+ spsrc <- as(st_read("spsrc.shp"), "Spatial")
77
+ sproads <- as(st_read("sproads.shp"), "Spatial")
74
78
 
75
79
 
76
80
  ###################################################
@@ -92,9 +96,9 @@
92
96
  # spatstat 1.45-0 enforce a tighter spacing of the r vector 2016-03-22
93
97
  # r <- seq(0, sqrt(2)/6, by = 0.005)
94
98
  r <- seq(0, sqrt(2)/6, by = 0.001)
95
- envjap <- envelope(as(spjpines1, "ppp"), fun=Gest, r=r, nrank=2, nsim=99)
96
- envred <- envelope(as(spred, "ppp"), fun=Gest, r=r, nrank=2, nsim=99)
97
- envcells <- envelope(as(spcells, "ppp"), fun=Gest, r=r, nrank=2, nsim=99)
99
+ envjap <- envelope(as.ppp(st_as_sf(spjpines1)), fun=Gest, r=r, nrank=2, nsim=99)
100
+ envred <- envelope(as.ppp(st_as_sf(spred)), fun=Gest, r=r, nrank=2, nsim=99)
101
+ envcells <- envelope(as.ppp(st_as_sf(spcells)), fun=Gest, r=r, nrank=2, nsim=99)
98
102
  Gresults <- rbind(envjap, envred, envcells)
99
103
  Gresults <- cbind(Gresults,
100
104
  y=rep(c("JAPANESE", "REDWOOD", "CELLS"), each=length(r)))
@@ -118,9 +122,9 @@
118
122
  ### code chunk number 24: sppa.Rnw:661-670 (eval = FALSE)
119
123
  ###################################################
120
124
  set.seed(30)
121
- Fenvjap<-envelope(as(spjpines1, "ppp"), fun=Fest, r=r, nrank=2, nsim=99)
122
- Fenvred<-envelope(as(spred, "ppp"), fun=Fest, r=r, nrank=2, nsim=99)
123
- Fenvcells<-envelope(as(spcells, "ppp"), fun=Fest, r=r, nrank=2, nsim=99)
125
+ Fenvjap<-envelope(as.ppp(st_as_sf(spjpines1)), fun=Fest, r=r, nrank=2, nsim=99)
126
+ Fenvred<-envelope(as.ppp(st_as_sf(spred)), fun=Fest, r=r, nrank=2, nsim=99)
127
+ Fenvcells<-envelope(as.ppp(st_as_sf(spcells)), fun=Fest, r=r, nrank=2, nsim=99)
124
128
  Fresults<-rbind(Fenvjap, Fenvred, Fenvcells)
125
129
  Fresults<-cbind(Fresults,
126
130
  y=rep(c("JAPANESE", "REDWOOD", "CELLS"), each=length(r)))
@@ -178,7 +182,8 @@
178
182
  bwq
179
183
 
180
184
  #Spatstat code
181
- mserw<-bw.diggle(as(spred, "ppp"))
185
+ spred_ppp <- as.ppp(st_as_sf(spred))
186
+ mserw<-bw.diggle(spred_ppp)
182
187
  bw<-as.numeric(mserw)
183
188
  bw
184
189
 
@@ -198,7 +203,8 @@
198
203
  ###################################################
199
204
  library(splancs)
200
205
  poly <- as.points(list(x = c(0, 0, 1, 1), y = c(0, 1, 1, 0)))
201
- sG <- Sobj_SpatialGrid(spred, maxDim=100)$SG
206
+ library(stars)
207
+ sG <- as(as(st_as_stars(st_bbox(st_as_sfc(spred)), n=100*100), "Spatial"), "SpatialGrid")
202
208
  grd <- slot(sG, "grid")
203
209
  summary(grd)
204
210
  k0 <- spkernel2d(spred, poly, h0=bw, grd)
@@ -215,13 +221,13 @@
215
221
  ###################################################
216
222
  cc <- coordinates(kernels)
217
223
  xy<-list(x=cc[,1], y=cc[,2])
218
- k4<-density(as(spred, "ppp"), .5*bw, dimyx=c(100, 100), xy=xy)
224
+ k4<-density(spred_ppp, .5*bw, dimyx=c(100, 100), xy=xy)
219
225
  kernels$k4<-as(k4, "SpatialGridDataFrame")$v
220
- k5<-density(as(spred, "ppp"), .5*.05, dimyx=c(100, 100), xy=xy)
226
+ k5<-density(spred_ppp, .5*.05, dimyx=c(100, 100), xy=xy)
221
227
  kernels$k5<-as(k5, "SpatialGridDataFrame")$v
222
- k6<-density(as(spred, "ppp"), .5*.1, dimyx=c(100, 100), xy=xy)
228
+ k6<-density(spred_ppp, .5*.1, dimyx=c(100, 100), xy=xy)
223
229
  kernels$k6<-as(k6, "SpatialGridDataFrame")$v
224
- k7<-density(as(spred, "ppp"), .5*.15, dimyx=c(100, 100), xy=xy)
230
+ k7<-density(spred_ppp, .5*.15, dimyx=c(100, 100), xy=xy)
225
231
  kernels$k7<-as(k7, "SpatialGridDataFrame")$v
226
232
  summary(kernels)
227
233
 
@@ -305,9 +311,9 @@
305
311
  ### code chunk number 43: sppa.Rnw:1420-1430
306
312
  ###################################################
307
313
  set.seed(30)
308
- Kenvjap<-envelope(as(spjpines1, "ppp"), fun=Kest, r=r, nrank=2, nsim=99)
309
- Kenvred<-envelope(as(spred, "ppp"), fun=Kest, r=r, nrank=2, nsim=99)
310
- Kenvcells<-envelope(as(spcells, "ppp"), fun=Kest, r=r, nrank=2, nsim=99)
314
+ Kenvjap<-envelope(as.ppp(st_as_sf(spjpines1)), fun=Kest, r=r, nrank=2, nsim=99)
315
+ Kenvred<-envelope(as.ppp(st_as_sf(spred)), fun=Kest, r=r, nrank=2, nsim=99)
316
+ Kenvcells<-envelope(as.ppp(st_as_sf(spcells)), fun=Kest, r=r, nrank=2, nsim=99)
311
317
  Kresults<-rbind(Kenvjap, Kenvred, Kenvcells)
312
318
  Kresults<-cbind(Kresults,
313
319
  y=rep(c("JAPANESE", "REDWOOD", "CELLS"), each=length(r)))
@@ -333,10 +339,10 @@
333
339
  ###################################################
334
340
  bwasthma<-.06
335
341
 
336
- pppasthma<-as(spasthma, "ppp")
337
- pppasthma$window<-as(spbdry, "owin")
342
+ pppasthma<-as.ppp(st_as_sf(spasthma["Asthma"]))
343
+ pppasthma$window<-as.owin(st_as_sf(spbdry))
338
344
 
339
- marks(pppasthma)<-relevel(pppasthma$marks$Asthma, "control")
345
+ marks(pppasthma)<-relevel(pppasthma$marks, "control")
340
346
 
341
347
 
342
348
 
@@ -362,9 +368,9 @@
362
368
  ### code chunk number 54: sppa.Rnw:1880-1888
363
369
  ###################################################
364
370
 
365
- spkratio0<-as(kcases, "SpatialGridDataFrame")
371
+ spkratio0<-image2Grid(list(x=kcases$xcol, y=kcases$yrow, z=t(kcases$v)))
366
372
  names(spkratio0)<-"kcases"
367
- spkratio0$kcontrols<-as(kcontrols, "SpatialGridDataFrame")$v
373
+ spkratio0$kcontrols<-image2Grid(list(x=kcontrols$xcol, y=kcontrols$yrow, z=t(kcontrols$v)))$z
368
374
  spkratio<-as(spkratio0, "SpatialPixelsDataFrame")
369
375
 
370
376
  spkratio$kratio <- spkratio$kcases/spkratio$kcontrols
@@ -393,7 +399,7 @@
393
399
  kcasesrel <- density(casesrel, bwasthma)
394
400
  kcontrolsrel <- density(controlsrel, bwasthma)
395
401
  kratiorel <- eval.im(kcasesrel/kcontrolsrel)
396
- rlabelratio[i,] <- as(as(kratiorel, "SpatialGridDataFrame"), "SpatialPixelsDataFrame")$v
402
+ rlabelratio[i,] <- as(image2Grid(list(x=kratiorel$xcol, y=kratiorel$yrow, z=t(kratiorel$v))), "SpatialPixelsDataFrame")$z
397
403
  pvaluemap <- pvaluemap + (spkratio$kratio < rlabelratio[i,])
398
404
  }
399
405
 
@@ -415,7 +421,10 @@
415
421
  spkratio$pvaluemap <- (pvaluemap+1)/(niter+1)
416
422
  imgpvalue <- as.image.SpatialGridDataFrame(spkratio["pvaluemap"])
417
423
  clpvalue <- contourLines(imgpvalue, levels=c(0,.05, .95, 1))
418
- cl <- ContourLines2SLDF(clpvalue)
424
+ xx = lapply(clpvalue, function(x) cbind(x$x, x$y)) # avoiding maptools::ContourLines2SLDF
425
+ xxx = split(xx, sapply(clpvalue, "[[", "level"))
426
+ x2a = st_sfc(lapply(xxx, st_multilinestring))
427
+ cl = as(st_as_sf(x2a, level=names(xxx)), "Spatial")
419
428
 
420
429
 
421
430
  ###################################################
@@ -461,8 +470,7 @@
461
470
  ### code chunk number 71: sppa.Rnw:2293-2295
462
471
  ###################################################
463
472
  rr<-relrisk(pppasthma, bwasthmap)
464
- spkratio$prob<-as(as(rr, "SpatialGridDataFrame"), "SpatialPixelsDataFrame")$v
465
-
473
+ spkratio$prob<-as(image2Grid(list(x=rr$xcol, y=rr$yrow, z=t(rr$v))), "SpatialPixelsDataFrame")$z
466
474
 
467
475
  ###################################################
468
476
  ### code chunk number 72: sppa.Rnw:2307-2317
@@ -727,4 +735,9 @@
727
735
  ## lines(r, envKIcov$lo-envKIcov$mmean, lty=2)
728
736
  ## lines(r, envKIcov$hi-envKIcov$mmean, lty=2)
729
737
 
738
+ (sI <- sessionInfo()) # check: no sp?
739
+
740
+ "rgdal" %in% c(names(sI$otherPkgs), names(sI$loadedOnly))
741
+ "rgeos" %in% c(names(sI$otherPkgs), names(sI$loadedOnly))
742
+ "maptools" %in% c(names(sI$otherPkgs), names(sI$loadedOnly))
730
743