Set working directory and read

setwd("/Users/annapeterson/Desktop/Classes/GEOG6000/Lab13")

#Read ozone data, change to WGS84, and set geometry
ozone = st_read("/Users/annapeterson/Desktop/Classes/GEOG6000/Lab13/ozone.shp",
                quiet = TRUE)
st_crs(ozone) = 4326
ozone_geom = st_geometry(ozone)

#Read in lake, state, and location data

lk = st_read("/Users/annapeterson/Desktop/Classes/GEOG6000/Lab13/ne_50m_lakes.shp",
             quiet = TRUE)
lk_geom = st_geometry(lk)
st = st_read("/Users/annapeterson/Desktop/Classes/GEOG6000/Lab13/ne_50m_admin_1_states_provinces_shp.shp",
             quiet = TRUE)
st_geom = st_geometry(st)
loc = st_read("/Users/annapeterson/Desktop/Classes/GEOG6000/Lab13/ne_50m_populated_places.shp",
              quiet = TRUE)

Plot Sample Variogram

We can change the width which is the distance interval the semivariance is calculated. When we increase the width we’re able to see some things that were likely missed because it was either set too low or too high. Too low would make you miss some things and too high would create too much smoothing.
Variogram Model Fitted

modNugget = 0.08
modRange = nrow(ozone)
modSill = 0.07

ozone_var = variogram(thres ~ 1, ozone)

ozone_vgm1 = vgm(psill = modSill, 
                "Sph", 
                range = modRange, 
                nugget = modNugget)

ozone_vgm2 = fit.variogram(ozone_var, ozone_vgm1)


Plotting Cir vs Sph is fairly subtle if you’re not paying attention, but I still think the Sph model fits better.

Map of Predicted Probabilities

mybbox = st_bbox(ozone)
pred_grid = st_as_stars(mybbox,
                        xlim = c(-94, -82),
                        ylim = c(36,45),
                        dx = 0.1,
                        dy = 0.1)
st_crs(pred_grid) = 4326
ozone_vgm2 = fit.variogram(ozone_var, ozone_vgm1)
ozone_ik = krige(thres ~ 1,
                 ozone,
                 pred_grid,
                 ozone_vgm2,
                 nmax = 40)
## [using ordinary kriging]


You can see pockets of concentrated ozone which I suspect are near cities. But to confirm, I should plot it again. I will use methods we learned from lab 10 with tm_shape.

That took forever to pick a palette for the ozone that would mildly pop on the continuous color! Looks like we have a really high concentration of Ozone around Milwaukee and Chicago. This makes sense because from the EPA’s website, “cars, power plants, industrial boilers, refineries, chemical plants, and other sources chemically react in the presence of sunlight.”