Mosquito Mapping

’’’’’‘’Mosquito Count Pre-processing’’’’’’’’’

This notebook aims to find “dubplicated” GPS points within a radius of 30m.

library(raster)
library(dplyr)
library(sf)
library(rgdal)

1. Load the raster layer as a reference to locate duplicated GPS points

# NDVI raster layer at 30m 
NDVI<-raster("F:\\Research\\Research_Cooperation\\ILRI_Mosquito_Mapping\\Mosquito\\Cover_Ratio\\Innitial\\NDVI_Mean_1km.tif")
NAvalue(NDVI)<-0
# Mosquito count shapefile data
Shapefile<-shapefile("F:\\Research\\Research_Cooperation\\ILRI_Mosquito_Mapping\\Mosquito\\Mosquito_Sampling\\WGS84_Points.shp")
# Convert Mosquito shapefile to sf object

Mos_sf<-st_as_sf(Shapefile)

row.names(Mos_sf) <- NULL # Reset index 

# Plot the shapefile and NDVI values
par(mfrow=c(1,2),pty = "s")
plot(NDVI, main="Mean NDVI Values")
plot(Shapefile, main="Mosquito Sampling Location", pch=19)

2. Extract raster position ID (row and column)

# Extract the raster column and row possition 
Pos_ID<-cellFromXY(NDVI, Shapefile)
# Convert Mosquito count to integer 

Mos_sf$Sum_Cx<-as.integer(Mos_sf$Sum_Cx)
# Merge the Pos_ID to the shapefile

Mos_sf$Ras_ID<-Pos_ID
# Convert Ras_ID to factor 
Mos_sf$Ras_ID<-as.factor(Mos_sf$Ras_ID)

3. Calculate the total of mosquito at each position

Mos_uni<- Mos_sf %>% group_by(Ras_ID)%>% summarise(Sum_total=sum(Sum_Cx))
# Reorder the raster ID column
Mos_uni$Ras_ID<-as.integer(as.character(Mos_uni$Ras_ID)) 

Mos_uni<-Mos_uni %>% arrange(Ras_ID)

# Drop duplicated raster ID in Mos_sf

Mos_drop<-Mos_sf[!duplicated(Mos_sf$Ras_ID),]

Mos_drop$Ras_ID<-as.integer(as.character(Mos_drop$Ras_ID)) 

Mos_final<- Mos_drop %>% arrange(Ras_ID)

4. Merge total count to shapefile

Mos_final$Counts<-Mos_uni$Sum_total

# Select only column of interest 
Mos_final<-Mos_final %>% select(ID, Lat,Long, District,Counts)

head(Mos_final)

Simple feature collection with 6 features and 5 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: 105.6406 ymin: 21.13722 xmax: 105.6461 ymax: 21.14083
CRS:            +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
        ID      Lat     Long   District Counts                  geometry
1 02 - 101 21.14083 105.6461 Dan Phuong     32 POINT (105.6461 21.14083)
2 02 - 098 21.13972 105.6450 Dan Phuong      4  POINT (105.645 21.13972)
3 02 - 100 21.13944 105.6458 Dan Phuong     43 POINT (105.6458 21.13944)
4 02 - 093 21.13806 105.6428 Dan Phuong      0 POINT (105.6428 21.13806)
5 02 - 096 21.13722 105.6406 Dan Phuong     11 POINT (105.6406 21.13722)
6 02 - 095 21.13722 105.6411 Dan Phuong      0 POINT (105.6411 21.13722)

• Another way to get the same output in Section 4.

Shape_final<-st_join(Mos_final,Mos_uni,by="Ras_ID",left=T)

although coordinates are longitude/latitude, st_intersects assumes that they are planar

head(Shape_final)

Simple feature collection with 6 features and 7 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: 105.6406 ymin: 21.13722 xmax: 105.6461 ymax: 21.14083
CRS:            +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
        ID      Lat     Long   District Counts  Ras_ID Sum_total                  geometry
1 02 - 101 21.14083 105.6461 Dan Phuong     32 2465648        32 POINT (105.6461 21.14083)
2 02 - 098 21.13972 105.6450 Dan Phuong      4 2476512         4  POINT (105.645 21.13972)
3 02 - 100 21.13944 105.6458 Dan Phuong     43 2479232        43 POINT (105.6458 21.13944)
4 02 - 093 21.13806 105.6428 Dan Phuong      0 2495523         0 POINT (105.6428 21.13806)
5 02 - 096 21.13722 105.6406 Dan Phuong     11 2503666        11 POINT (105.6406 21.13722)
6 02 - 095 21.13722 105.6411 Dan Phuong      0 2503668         0 POINT (105.6411 21.13722)

5. Write the shapefile to our local computer

# Set the path to save
path_save<-"F:\\Research\\Research_Cooperation\\ILRI_Mosquito_Mapping\\Mosquito\\Mosquito_Sampling"

# st_write(Mos_final,dsn = path_save,"Mosquito_Final.shp",driver="ESRI Shapefile",delete_layer=T)