library(terra)
library(sf)
library(stars)
library(ggplot2)
library(viridis)
#Step1: Read the raster
b <- terra::rast("./temp/cru_ts4.07.1901.2022.tmp.dat.nc")
#Step2: Extract time information
time_values <- time(b)
#Step3: Turning strings into time
date_objects <- as.Date(time_values)
# Step 4: Select only the dates from the year 1972 and 2022
dates_1972 <- date_objects[format(date_objects, "%Y") == "1972"]
dates_2022 <- date_objects[format(date_objects, "%Y") == "2022"]
# Step 5: Find the index of the selected dates in the time_values
selected_indices_1972 <- which(time_values %in% dates_1972)
selected_indices_2022 <- which(time_values %in% dates_2022)
# Step 6: Extract raster values for the selected dates
selected_rasters_1972 <- b[[selected_indices_1972]]
selected_rasters_2022 <- b[[selected_indices_2022]]
# Step 7: Extract layer names
layer_names_1972 <- names(selected_rasters_1972)
layer_names_2022 <- names(selected_rasters_2022)
# Step 8: Find indices of layers that start with "tmp"
indices_tmp_1972 <- grep("^tmp", layer_names_1972)
tmp_layers_1972 <- selected_rasters_1972[[indices_tmp_1972]]
indices_tmp_2022 <- grep("^tmp", layer_names_2022)
tmp_layers_2022 <- selected_rasters_2022[[indices_tmp_2022]]
# Step 9: Changing raster names to temperature
# Ensure that the number of names matches the number of layers
if(length(indices_tmp_1972) == length(dates_1972)) {
names(tmp_layers_1972) <- format(dates_1972, "%Y-%m-%d")
}
if(length(indices_tmp_2022) == length(dates_2022)) {
names(tmp_layers_2022) <- format(dates_2022, "%Y-%m-%d")
}
#Step10: Turning the spatraster into a list of raster
stars_list_1972 <- lapply(tmp_layers_1972, function(raster_obj) st_as_stars(raster::raster(raster_obj)))
stars_list_2022 <- lapply(tmp_layers_2022, function(raster_obj) st_as_stars(raster::raster(raster_obj)))
#Step11: Use lapply to apply st_apply to each stars object to calculate the mean
mean_list_1972 <- lapply(stars_list_1972, function(stars_obj) st_apply(stars_obj, c("x", "y"), mean))
mean_list_2022 <- lapply(stars_list_2022, function(stars_obj) st_apply(stars_obj, c("x", "y"), mean))
#Step12: Combine the list of results into a single stars object
mean_stars_1972 <- do.call(stars::st_as_stars, mean_list_1972)
mean_stars_2022 <- do.call(stars::st_as_stars, mean_list_2022)
#Step13:Plotting Global Average Temperature in 1972
fig <-ggplot() +
geom_stars(data = mean_stars_1972) +
scale_fill_viridis_c( na.value = NA) +
theme_bw() +
ggtitle("1972") +
labs(fill = "" )+
theme(
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank())
#Step14:Plotting Global Average Temperature in 2022
fig2<-ggplot()+
geom_stars(data=mean_stars_2022)+
scale_fill_viridis_b(na.value=NA)+
theme_bw()+
ggtitle(2022)+
labs(fill = "" )+
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank())Rising Temperatures and Droughts
Focus on New Zealand
1 Introduction
Research Question: How has Climate Change Affected Droughts in New Zealand?
This analysis will follow a comparative approach, focusing on the years 1972 and 2022. The choice of these two years is strategic; 1972 marks a landmark year in environmental politics, while 2022 provides a year for which sufficient data were found, yet is not too dated.
First the year 1972 marks the UN Conference on the Human Environment in Stockholm. This was the first world conference to prioritize the environment as a major issue. Although the main focus of the conference was not climate change, it set the stage for future discussions on the topic.
Overtime, Many other environmental conventions are worth mentioning, such as the 1992 Rio Summit and the 2000 Millennium Summit, which eventually led to the publication of the Millennium Development Goals, to be achieved by 2015.
Thus, 2015 is a crucial year for environmental politics since two main milestones of environmental politics were also achieved:
- Paris Agreement: Adopted on December 12, 2015, at COP21 in Paris, it is landmark international treaty within the UNFCCC to limit global warming to well below 2 degrees Celsius, preferably to 1.5 degrees Celsius, compared to pre-industrial levels
- Sustainable Development Goals: 17 global goals adopted by UN Member States, aimed to be achieved by 2030.
For the purpose of this presentation, what we should have in mind is Goal 13. 
2 Data
For this research, I have used a dataset of average temperatures provided in class, selecting only the years relevant to my study. Additionally, I employed New Zealand shapefiles from the GADM site at levels 0, 1, and 2 to accurately define geographical boundaries.
To measure drought conditions up to the year 2022, I utilized the Standardised Precipitation-Evapotranspiration Index (SPEI), which integrates both temperature and precipitation data, using data called Data to 2022
3 Global Context
Before delving into the core subject of the analysis, Average Temperature Change from 1972 to 2022 will be examined using world maps to display geospatial data,generated utilizing ‘stars’ for advanced raster handling and visualization.
3.1 Global Average Temperature in 1972 vs 2022
There has been a clear shift toward higher temperatures, particularly in the northern hemisphere. This includes a notable warming in the Arctic region, where temperatures have transitioned from cooler blues in 1972 to warmer greens and yellows by 2022, indicating significant overall warming.
3.2 Time Trend of Temperatures
This is further supported examining a the time trend of temperatures where blue bars represent years with cooler average temperatures and red bars indicate years with higher average temperatures.Earlier years exhibit a clustering of blue bars, denoting cooler temperatures.Over time, there is a noticeable increase in red bars, suggesting a warming trend from 1972 to 2022.
3.3 Context Relevance
Considering the analysis made thus far, acknowledging that temperature is increasing globally is vital since, according to a 2024 report of the European Commission called Consequences of Climate Change, higher temperatures increase the evaporation of water, which, along with the lack of precipitation, heightens the risk of severe droughts. Consequently, this leads us to the next step of the analysis.
4 Focus on New Zealand
This section will focus on the average temperature change, conducting a GIS analysis to determine whether the European Commission’s statement is applicable to New Zealand. This country has been selected for this analysis due to its historical susceptibility to drought conditions. Its unique geography and climate patterns have led to recurrent periods of water scarcity, making it an important case study for understanding the impact of rising temperatures on drought intensity.
The first step is to analyze whether the average temperature increased using zonal statistics for enhanced visualization, utilizing shapefiles of levels 0 and 2 alongside rasters to plot the average temperature. the viridis scale ‘B’ is employed here, ensuring a clear visual distinction in the temperature differences.
The 1972 map mainly features purple and dark blue shades, indicating cooler temperatures. The 2022 map shows a significant shift to green and yellow tones, especially in the northern regions, suggesting warmer average temperatures compared to 1972. Overall, there is an increase in average temperatures over the 50-year span between these two years.
Raster Difference
This raster difference map provides a detailed view of temperature change distribution. The southern South Island is shown in dark red, marking the highest temperature increase of up to 2.25 degrees Celsius over 50 years. Conversely, the northern South Island and much of the North Island, depicted in green to yellow, experienced smaller increases, ranging from approximately 1.25 to 1.75 degrees Celsius.
5 Regional Analysis
First, to simplify visualization and analysis, I plotted the centroids of New Zealand’s 16 regions. This allows the subsequent analysis to focus on regional data rather than district-level details.
New Zealand Regions 
This section will analyze ‘intensity’ as a drought metric, which is defined as severity per duration. Additionally, maps displaying other metrics such as frequency (number of events), duration (months), and severity (dryness level) can be produced upon request.
5.1 Droughts Intensity Compared
In 1972, we can observe a scattering of mostly light yellow and orange circles, suggesting lower to moderate drought intensities across the regions. Contrastingly, the 2022 map shows an increase in the number of orange circles and the emergence of dark red areas, especially in the northern part of the South Island and parts of the North Island, indicating a general increase in drought intensity over the 50-year period. This suggests that New Zealand has experienced a notable rise in drought intensity in duration, aligning with the global trend of increasing temperatures and climate change effects.
5.2 Mean Intensity Compared
Also, to assess the average drought intensity per region, I grouped the data by region name and calculated the mean intensity, excluding missing values for both 1972 and 2022.
The 1972 drought intensity map for New Zealand shows varying conditions, with the South Island displaying a gradient from light yellow in the north to blue in the south, indicating mild to moderate drought. In contrast, the North Island predominantly exhibits purple hues, signaling severe drought in some regions. The 2022 map reveals a sharper contrast, with extensive data gaps and concentrated areas of dark purple and blue in the southern South Island and parts of the North Island, denoting more intense drought conditions in these locations.
6 Conclusion
6.1 Summary
Taking into account the data presented over a 50-year period (1972-2022), temperature has indeed increased, despite efforts to combat climate change through emissions reduction. This upward temperature trend is evident through visualizations and time trend analyses. Consequently, it has been observed that temperatures in New Zealand have also risen. This is particularly significant given that New Zealand has historically been susceptible to droughts, which are now intensifying.
6.2 Solutions
This calls for an approach that would address the increasing high temperatures and increased levels of droughts in New Zealand. Thus, more appropriate water conservation measures have to be set up, which include harvesting rainwater and growing drought-tolerant crops. Besides, advanced irrigation technologies could enhance improved water use and decrease wastage.
Further, reforestation and rehabilitation of the wetlands will contribute very much to the capacity of the landscape to retain water. Furthermore, in order to reduce greenhouse gas emissions and moderate the impacts of climate change, appropriate policies should be applied to make sure that there are proper sustainable farming and use of renewable energy since these contributions will reflect upon them.
Finally, the resilience to the altering climate challenge is something that, if brought out, will have been due to shared efforts between the community and government support.
6.3 Policy Implications
As temperatures continue rising with drier conditions being recorded in the country, policies within New Zealand will have to shift further to resilient and sustainable paths. In the process, a better national water strategy has to be established, focusing on issues of conservation and efficiency. It would require making compulsory the installation of water-conserving equipment and fixtures in new building constructions.
Incentives would also be made for those practices in agriculture that conserve water and save soil moisture. In addition, some other measures to provide more financial support to the development of drought-resistant plant variety and sustainable farming technologies should be undertaken. It, therefore, requires an adjustment of the regulations in urban planning tp create water-permeable surfaces to prevent excessive runoff and contribute to topping up groundwater.
Finally, government investment in renewable sources of energy can enable us to reduce our carbon footprint.
7 Methodology
This section includes a description of the library and operations used in each section to create the maps showed in this document.
7.1 Global Context
To create the maps of Global Average Temperature between 1972 vs. 2022, spatial data manipulation and visualization operations were carried out.
Libraries used: - terra: Used for handling raster data effectively, this is a kind of spatial operations on the global temperature datasets.
sf: An array of tools was presented to handle and manipulate geometric and geographic data.
stars: This has helped in data handling where raster data has spatiotemporal arrays, aligning, and comparing data across different time frames.
Key functions:
- terra::rast: Generated raster objects from the temperature datasets which provide a basis in the conduction of spatial analyses.
- lapply: This function will apply certain functions over the list of datasets to streamline the data processing for both years.
- stars::st_as_stars: Converts the raster data into stars objects to treat and visualize temporal dimensions of the datasets easily.
7.2 Focus on New Zealand
To produce the focus on new Zealand map and the Raster Difference map,the following R libraries and functions were used to show the average temperature and difference in raster data:
Libraries Used:
- units: Handle and convert between different measurement units in spatial datasets for use in spatial operations.
Key functions:
- st_crop: This function is in the sf package and selects New Zealand by excluding other spatial data from the global context, hence it is important in regional analysis.
- raster::aggregate: To reduce resolution in the raster data for computation, this helps in dealing with larger datasets while keeping the desired details during computation.
- New Zealand shapefiles of level 0 and 2: These tools were important to depict the climatic conditions accurately over different regions in New Zealand.
7.3 Regional Analysis
In this section that looks at New Zealand Regions, the intensity and mean intensity were looked at using specific R libraries and functions to anlayze the spatial data distributions with more meaning in precision.
Libraries Used:
- dplyr: This package facilitated manipulation and summarization of the data, essential while dealing with large datasets, including preparing them for spatial analysis.
Key Functions:
st_area, set_units: This set of functions was used to calculate the area of geographic features within the provided New Zealand shapefiles and set the units of measurement that will establish valid area calculations from the dataset.
st_centroid: Used in calculating the centroid of shapes important in finding the average mean point of intensity in given regions.
Level 1 New Zealand shapefiles: These form an integral part of the detailed geographic boundary shapefile necessary for the regional breakdown to permit the carrying out of detailed spatial analysis.
8 Appendix
This section will include the code necessary to create all the maps and visual data in helpful in each section
8.1 Global Context
Global Average Temperature in 1972 vs 2022
Time Trend of Temperatures
# Reading the file
average_temp <- read.csv("./average_temp_1901_2022.csv")
# Subsetting to between 1972 and 2022 (inclusive)
average_temp_1972_2022 <- subset(average_temp, year >= 1972 & year <= 2022)
# Calculating the 1972-2022 average
average_temp_1972_2022_avg <- mean(average_temp_1972_2022$average_temperature)
# Calculate the deviations from the 1972-2022 average only for the years 1972 to 2022
average_temp_1972_2022$dev_1972_2022_avg <- average_temp_1972_2022$average_temperature - average_temp_1972_2022_avg
# Create a new column for bar color based on the condition, but only for the years 1972 to 2022
average_temp_1972_2022$bar_color <- ifelse(average_temp_1972_2022$dev_1972_2022_avg > 0, "red", "blue")
# Add the deviation and bar color columns to the original dataframe for years 1972 to 2022
average_temp <- merge(average_temp, average_temp_1972_2022[, c("year", "dev_1972_2022_avg", "bar_color")], by = "year", all.x = TRUE)
# Create the bar plot
v<-ggplot(average_temp, aes(x = year, y = dev_1972_2022_avg, fill = bar_color)) +
geom_bar(stat = "identity", color = "black") +
scale_fill_manual(values = c("blue", "red"), guide="none") +
theme_bw() +
labs(title = "Time Trend of Average Temperature",
x = "Year",
y = "Difference from 1972-2022 avearge (degrees C)")8.2 Focus on New Zealand
Focus on New Zealand
library(sf)
library(ggplot2)
library(units)
library(ggpubr)
# Reading borders at level 0 and 2
nz0 <- st_read(dsn="./gadm41_NZL_shp/gadm41_NZL_0.shp", quiet = TRUE)
nz2 <- st_read(dsn="./gadm41_NZL_shp/gadm41_NZL_2.shp", quiet = TRUE)
#Simplifying the shapefile for easier plotting
nz0<-st_simplify(nz0, dTolerance=0.02)
nz2<-st_simplify(nz2, dTolerance=0.005)
#Extracting only the relevant variable
nz2b<-subset(nz2, select = c("NAME_2"))
#Calculating the area
nz2b$area<-st_area(nz2b)
nz2b$area_sqkm<-set_units(nz2b$area, 'km^2')
nz2b$area_sqkm_numeric<-as.numeric(nz2b$area_sqkm)
#setting extents
min_lon_x <- 166
max_lon_x <- 179
min_lat_y <- - 48
max_lat_y <- -34
#1972
#Ensuring that the CRS is the same
st_crs(mean_stars_1972)<-st_crs(nz0)
# Crop the resampled raster data to NZ boundaries
cropped_1972<- st_crop(mean_stars_1972, nz0)
# Aggregating lights by polygon
agg_1972<- raster::aggregate(cropped_1972, nz2, mean, na.rm=TRUE)
f_1972<- ggplot()+
geom_stars(data=agg_1972)+
coord_sf(
xlim = c(min_lon_x , max_lon_x),
ylim = c(min_lat_y , max_lat_y)
) +
ggtitle("1972")+
labs(fill = "" ) +
scale_fill_viridis_b(na.value = NA)
#2022
# Ensuring that the CRS is the same
st_crs(mean_stars_2022)<-st_crs(nz0)
# Crop the resampled raster data to China's boundaries
cropped_2022<- st_crop(mean_stars_2022, nz0)
#Step2: Aggregating lights by polygon
agg_2022<- raster::aggregate(cropped_2022, nz2, mean, na.rm=TRUE)
f_2022<-ggplot()+
geom_stars(data=agg_2022)+
coord_sf(
xlim = c(min_lon_x , max_lon_x),
ylim = c(min_lat_y , max_lat_y)
) +
ggtitle("2022")+
labs(fill = "" ) +
scale_fill_viridis_b(na.value = NA)
# Comparing the Two Plots
mean_temp<-ggarrange(f_1972, f_2022, common.legend = TRUE)
arranged_mean <- annotate_figure(mean_temp,top = text_grob("Average Temperature 1972 vs.2022",
size = 14,
face = "bold"))Raster Difference
library(ggplot2)
library(viridis)
# Calculate the difference
diff_raster <- agg_2022 - agg_1972
# Plot the difference using ggplot and viridis theme
diff_map <- ggplot() +
geom_stars(data = diff_raster) +
scale_fill_viridis_c(option="turbo", name = "Temp.
2022-1972", na.value = NA) +
theme_bw()+
coord_sf(
xlim = c(min_lon_x , max_lon_x),
ylim = c(min_lat_y , max_lat_y)
) +
theme(legend.position = c(0, 1),
legend.justification = c(0, 1),
legend.background = element_rect(color = "black", size = 0.2)) 8.3 Regional Analysis
New Zealand Regions
library(sf)
sf_use_s2(FALSE)
library(ggplot2)
library(ggrepel)
# Loading New Zealand Shapefile
nz_cntry1 <- st_read(dsn="./gadm41_NZL_shp/gadm41_NZL_1.shp", quiet = TRUE)
# Simplify lines
nz1<-st_simplify(nz_cntry1, dTolerance = 0.05)
#Extracting only the relevant variable
nz1b<-subset(nz1, select = c("NAME_1"))
#Calculating the area
nz1b$area<-st_area(nz1b)
nz1b$area_sqkm<-set_units(nz1b$area, 'km^2')
nz1b$area_sqkm_numeric<-as.numeric(nz1b$area_sqkm)
library(ggrepel)
# Extract centroids
centroids <- st_centroid(nz1b)
# Convert to a simple data frame for ggplot
centroids_df <- as.data.frame(st_coordinates(centroids))
centroids_df$NAME_1 <- nz1b$NAME_1
centroids_df$area_sqkm_numeric <- nz1b$area_sqkm_numeric
k<-ggplot(nz1b, aes(fill = area_sqkm_numeric)) +
geom_sf() +
geom_label_repel(data = centroids_df, aes(x = X, y = Y, label = NAME_1),
colour = "black", size = 2, fill = "white",
box.padding = unit(0.6, "lines")) + # Adjust box.padding as necessary
coord_sf(xlim = c(min_lon_x, max_lon_x), ylim = c(min_lat_y, max_lat_y)) +
scale_fill_viridis_c(option = "viridis")+
theme_void()Droughts Intensity Compared
library(dplyr)
library(sf)
sf_use_s2(FALSE)
library(ggplot2)
#Step1: Reading droughts
droughts <- read.csv('./drought.data/drought-spei-1962-2022.csv')
#Step2: Turning the CSV file into an SF object
droughts_sf<-st_as_sf(droughts, coords = c("lon", "lat"))
#Step4: Associating the crs information for the NZ shapefile to points
droughts_sf <- st_set_crs(droughts_sf, st_crs(nz1))
#1972
droughts_1972 <- subset(droughts_sf, year == "1972")
#2022
droughts_2022 <- subset(droughts_sf, year == "2022")
# 1972: Plotting Intensity
plotting_nz<-ggplot() +
geom_sf(data = nz1b, fill = "white", color = "black", size = NA) + # Base map
geom_sf(data = droughts_1972, aes(size = intensity, color = intensity), alpha = 0.7) + # Drought data
scale_color_gradient2(
low = "blue", mid = "yellow", high = "red",
breaks = c(0.25, 0.50, 0.75),
name = "Intensity"
) +
scale_size_continuous(
breaks = c(0.25, 0.50, 0.75), # Size breaks matching color breaks
range = c(0.2, 3),
name = "Intensity"
) +
coord_sf(
xlim = c(min_lon_x, max_lon_x),
ylim = c(min_lat_y, max_lat_y)
) +
theme_minimal() +
labs(x = "Longitude", y = "Latitude", title = "1972: Drought Intensity in New Zealand") +
guides(
color = guide_legend(),
size = guide_legend()
)
# 2022: Plotting Intensity
plotting_22 <- ggplot() +
geom_sf(data = nz1b, fill = "white", color = "black", size = NA) + # Base map
geom_sf(data = droughts_2022, aes(size = intensity, color = intensity), alpha = 0.7) + # Drought data
scale_color_gradient2(
low = "blue", mid = "yellow", high = "red", # Adjusted midpoint and colors
breaks = c(0.25, 0.50, 0.75),
name = "Intensity"
) +
scale_size_continuous(
breaks = c(0.25, 0.50, 0.75), # Size breaks matching color breaks
range = c(0.2, 3),
name = "Intensity"
) +
coord_sf(
xlim = c(min_lon_x, max_lon_x),
ylim = c(min_lat_y, max_lat_y)
) +
theme_minimal() +
labs(x = "Longitude", y = "Latitude", title = "2022: Drought Intensity in New Zealand") +
guides(
color = guide_legend(),
size = guide_legend()
)Mean Intensity Compared
#1972
x1<-st_join(droughts_1972,nz1b)
xdf2<- st_drop_geometry(x1)
sdf3<-xdf2%>%
group_by(NAME_1)%>%
summarize(mean_int = mean(intensity, na.rm=TRUE))
fin_df72<-left_join(nz1b, sdf3, by = c('NAME_1'='NAME_1'))
#Plotting
fig3<- ggplot() +
geom_sf(data = fin_df72, aes(fill = mean_int), linewidth = 0.07) +
scale_fill_viridis_b (na.value = NA, name = "") +
labs(title = "1972", fill = "Drought Intensity") +
coord_sf(
xlim = c(min_lon_x , max_lon_x),
ylim = c(min_lat_y , max_lat_y)
) +
labs(x = "", y = "")+
theme_bw()
#2022
x<-st_join(droughts_2022,nz1b)
xdf<- st_drop_geometry(x)
sdf2<-xdf%>%
group_by(NAME_1)%>%
summarize(mean_int_2 = mean(intensity, na.rm=TRUE))
fin_df<-left_join(nz1b, sdf2, by = c('NAME_1'='NAME_1'))
#Plotting
fig4<-ggplot() +
geom_sf(data = fin_df, aes(fill = mean_int_2), linewidth = 0.07) +
scale_fill_viridis_b(name = "", na.value = NA) +
labs(title = "2022", fill = "Drought Intensity") +
coord_sf(
xlim = c(min_lon_x , max_lon_x),
ylim = c(min_lat_y , max_lat_y)
) +
labs(x = "", y = "")+
theme_bw()