Extract () augmente la somme des valeurs?

Question

J'utilise la grille Landscan pour calculer le nombre de personnes dans les communes vietnamiennes. J'ai d'abord masqué les données avec la limite GADM. Ensuite, j'ai utilisé extract () pour calculer la population par commune et enregistré sous forme de fichier sp. J'ai converti le fichier sp en une trame de données. Le code est ci-dessous.

pop <- raster("w001001.adf") 
crop <- crop(pop, gadm)
mask1 <- mask(crop, gadm) # this is so that the next step is a bit quicker
vn <- extract(mask1, gadm, fun=sum, na.rm=T, sp=T)
dta <- as.data.frame(vn) 

sum(dta$mask1) # sum of data frame ~ 82 millions
cellStats(mask1, sum) # sum of raster value ~ 95 millions

Je pense que la somme du raster ou de la base de données devrait me donner la population entière du Vietnam et donc être identique. Cependant, les deux sommes sont différentes. La somme des valeurs raster me donne 82 millions (ce qui est probablement correct) et la somme de la trame de données donne 95 millions de personnes ... Est-ce que quelqu'un sait ce que je fais mal?

Mettre à jour:

cell <-extract(mask1, gadm, cellnumbers=T, na.rm=T)
> sum(duplicated(cell))
[1] 7
> cell[duplicated(cell)]
[[1]]
    cell value
 1332229 47417

[[2]]
    cell value
 1332228 38578

[[3]]
    cell value
 1332227 41256

[[4]]
    cell value
 1329593 27005

[[5]]
    cell value
 1332233 41888

[[6]]
       cell value
[1,] 248323 29696
[2,] 249202 28011

[[7]]
   cell value
 248324 21301

Mise à jour 2:

Les données Landscan sont disponibles sur https://landscan.ornl.gov/landscan-datasets. Cependant, il existe un processus d'approbation qui prend quelques jours avant que les données puissent être téléchargées. J'essaie de faire exactement la même chose avec la valeur F162004 DMSP Nightlight qui peut être téléchargée à partir dehttps://ngdc.noaa.gov/eog/dmsp/downloadV4composites.html.

Le fichier de limites GADM peut être chargé directement dans r

library(raster)
gadm <- getData('GADM', country='VN', level=3) # VN

Answer 1

C'est parce que vous avez beaucoup de petits polygones et que, par défaut, il small=TRUEest défini. Dans ce cas, lorsqu'un polygone est si petit (ou peut-être de forme étrange) qu'il ne couvre pas le centre d'un carré de la grille, il extractfait autre chose (ce n'est pas clair exactement quoi) pour obtenir une valeur pour ce polygone.

Avec les données worldpop que j'obtiens:

> vn <- extract(mask1, gadm, fun=sum, na.rm=T, sp=T)
> sum(pop[],na.rm=TRUE)
[1] 99488591
> sum(vn$vnm_ppp_2020_1km_Aggregated)
[1] 118414408

qui est comme ce que vous avez, mais utiliser small=FALSEfait ceci:

> vn_nosmall <- extract(mask1, gadm, fun=sum, na.rm=T, sp=T, small=FALSE)
> sum(vn_nosmall$vnm_ppp_2020_1km_Aggregated)
[1] NA
> sum(vn_nosmall$vnm_ppp_2020_1km_Aggregated,na.rm=TRUE)
[1] 98847054

Premièrement, les petits polygones ont NApour valeur. Je peux savoir quelle est leur population par défaut, et cela totalise jusqu'à:

> sum(vn$vnm_ppp_2020_1km_Aggregated[
      is.na(vn_nosmall$vnm_ppp_2020_1km_Aggregated)
    ])
[1] 19567354

ce qui correspond à peu près à la différence entre le total raster et le total le plus élevé. L'excédent est un double comptage à partir de ces petits polygones.

help(extract)a quelques options qui peuvent résoudre ce problème, voir weightspar exemple, mais cela implique l'intersection de chaque polygone avec la grille pour additionner tous les petits bits et peut être très lent.

Vous pouvez également créer votre grille à une échelle plus fine afin que les petits polygones couvrent au moins un centre carré de grille. Vous pouvez le faire en divisant chaque cellule en quatre (par exemple en utilisant disaggregate), mais n'oubliez pas de diviser toutes les variables de comptage par quatre (mais pas les variables de densité qui sont constantes par unité de surface).