GY667 Assignment 2 - Ocean Indices, SST and SLP

Introduction

As a whole, ocean circulation is result of the contribution of various atmosphere-ocean system interactions which produce both locally specific and globally detectable influences on the general distribution of heat from the equator to the poles (NOAA, 2019; Nelson, 2016; Munk, 1950). Most notable, deep water formation and meridional overturning circulation. In considering changes in climate indices as sea surface temperature (SSTs) that are often quoted among the resulting impacts of anthropogenic climate change (Pachauri and Meyer, 2014:7), it is necessary to consider these other earth system dynamics , and in particular, oscillations which may impact variability in these indices not just over time, but over space too (Peterson and White, 1998; Bader and Latif, 2003).

In this assignment, SST data sourced from the Hadley Centre Sea Ice and Sea Surface Temperature data set (HadISST) will be utilised in order to assess potential spatially varying correlations in three ocean oscillations/climate indices:

1. Atlantic Multidecadal Oscillation/Variability
2. Pacific Decadal Oscillation
3. El Nino/Southern Oscillation

Regression analysis will also be carried out for each indice to identify patterns. Following the documentation of the methods used to produce each of these correlation and regression maps, a discussion of the emergent patterns will follow in light of the above considerations regarding ocean circulation.

Pre-requisite tasks

Setting the working directory:

setwd("F:/Assignment 2 rev01")

Ensuring R Studio is cleared prior to analysis:

graphics.off()
rm(list=ls())

Linking the folder containing the nc files.

path <- file.path(getwd(),"Data")

Loading the required libraries:

library(sp)

## Warning: package 'sp' was built under R version 3.5.3

library(raster) 
library(ncdf4)
library(rgdal)

## rgdal: version: 1.3-4, (SVN revision 766)
##  Geospatial Data Abstraction Library extensions to R successfully loaded
##  Loaded GDAL runtime: GDAL 2.2.3, released 2017/11/20
##  Path to GDAL shared files: C:/Users/Jordan/Documents/R/win-library/3.5/rgdal/gdal
##  GDAL binary built with GEOS: TRUE 
##  Loaded PROJ.4 runtime: Rel. 4.9.3, 15 August 2016, [PJ_VERSION: 493]
##  Path to PROJ.4 shared files: C:/Users/Jordan/Documents/R/win-library/3.5/rgdal/proj
##  Linking to sp version: 1.3-1

library(RColorBrewer)
library(lattice)
library(latticeExtra)
library(reshape2)
library(maps)

## Warning: package 'maps' was built under R version 3.5.3

In order to map the results of the global correlation and regression patterns, a shapefile of the coastline features of the globe was downloaded from the Natural Earth Data website, and projected to WGS 84:

world.coast <- readOGR(dsn = file.path(path,"ne_110m_coastline"), layer = "ne_110m_coastline")

## OGR data source with driver: ESRI Shapefile 
## Source: "F:\Assignment 2 rev01\Data\ne_110m_coastline", layer: "ne_110m_coastline"
## with 134 features
## It has 3 fields
## Integer64 fields read as strings:  scalerank

world.coast <- spTransform(world.coast, CRS("+proj=longlat +datum=WGS84"))

Methods & Results

Assignment question: ‘Download the timeseries’

Below, each of the climate indice as well as the sea level pressure data were downloaded and imported. From the KNMI Climate Explorer website, the HADISST versions of each indice were chosen in order to compare the correlations within the bounds of the HADISST methodology. This was deemed appropriate due to the potential for other datasets potentially having different errors or biases which could also affect correlation results.

HAD <-nc_open(file.path(path,'HadISST_sst.nc')) #Workshop 3 #1870 starts?
AMO <-nc_open(file.path(path, 'iamo_hadsst_ts.nc')) #AMOHADSST #1850-now
PDO <-nc_open(file.path(path, 'ipdo_hadsst3.nc')) #PDOHADSST #1850-now
ENSO <-nc_open(file.path(path, 'ihadisst1_nino4a.nc')) #HADISST #1870-now, HadISST1
SLP <-nc_open(file.path(path, 'slp.mon.mean.nc')) #From Moodle

Assignment question: ‘Calculate the correlation and regression patterns associated with each for SST and SLP’.

Hadley Centre Sea Ice and Sea Surface Temperature data set (HadISST)

In order to complete the this (and all four) objectives of this assignment, the HADISST data would need to be first imported and formatted accordingly. The proceeding code illustrates the process of ensuring that annually averaged SSTs could be used in the correlation calculations with each indices’.

First, by using print, it was identified that the SST data started in 1870, and was given in days since then. The next steps involved accessing the time variable in the nc data-set, which was also displayed in the print results.

H.time<-ncvar_get(HAD, 'time') 

In order to extract the years from which corresponding SST could be analysed, the ‘year’ sub-string was identified in the character string and extracted:

H.tunits<-ncatt_get(HAD,"time",attname="units") 
H.tustr<-strsplit(H.tunits$value, " ") 
H.date<-as.character(as.Date(H.time,origin=unlist(H.tustr)[3]))

The nc file also contained latitude and longitude, as well as the SST data. All of these variables were accessed:

lat<-ncvar_get(HAD, 'latitude')
lon<-ncvar_get(HAD, 'longitude')

H.sst <- ncvar_get(HAD, 'sst')

As the data-set contained metadata on the coding of missing values as well as time, the following code was used to account for this and to ensure it doesn’t affect the creation of annual averages of SST globally:

H.fillvalue <- ncatt_get(HAD,"sst","_FillValue")  
H.sst[H.sst==H.fillvalue$value] <- NA 
H.missvalue <- ncatt_get(HAD,"sst","missing_value")
H.sst[H.sst==H.missvalue$value] <- NA 
H.sst[H.sst==-1000] <- NA

By isolating the ‘year’ variable, the ‘colMeans’ and ‘aggregate’ functions were used to create a global mean for the SST variable and then to group these together for their associated year, respectively:

H.year <- format(as.Date(H.date, format="%Y-%m-%d"),"%Y")
H.gmean <- colMeans(H.sst, na.rm = TRUE, dims=2)
H.annmean <- aggregate(H.gmean,by=list(H.year),FUN=mean,na.rm=TRUE)

However, to identify average sea surface temperature, the rowMeans function was used.

H.avsst = rowMeans(H.sst,na.rm=FALSE,dims=2)

With the variables of interest now calculated, the RColorBrewer package could be used to set up a colour palette for mapping the SSTs on a levelplot as a vector:

colors <- rev(brewer.pal(10, "RdYlBu"))
pal <- colorRampPalette(colors)

grid <- expand.grid(x=lon, y=lat)
grid$H.avsst <- as.vector(H.avsst)

levelplot(H.avsst~x*y,grid,col.regions = pal(100),
          xlab='Longitude',ylab='Latitude',main='Average SST') + 
  layer(sp.lines(world.coast))

However, in order to make annual averages to account for long-term variability, a loop was used which accounted for the length of the HADISST data and infilled annual means:

yrs <- H.annmean$Group.1 
nyr <- length(yrs)

H.asst <- array(NA,c(dim(lon),dim(lat),nyr))

for (k in 1:nyr) {
  H.asst[,,k] <- rowMeans(H.sst[,,H.year==yrs[k]],na.rm=FALSE,dims=2)
}

Creating a new column within the grid dataframe would then allow for the mapping of the annually averaged data:

grid$H.an_avsst <- as.vector(rowMeans(H.asst,na.rm=FALSE,dims=2))

SSTMap <- levelplot(H.an_avsst~x*y, data=grid,col.regions = pal(100),
          xlab='Longitude',ylab='Latitude',main='Annually Averaged SST') + 
  layer(sp.lines(world.coast))

SSTMap

Lastly, to assess variability, the global mean was removed from each of these years:

H.gmean <- colMeans(H.asst, na.rm = TRUE, dims=2)
for (k in 1:nyr){
  H.asst[,,k]<-H.asst[,,k]-matrix(H.gmean[k],length(lon),length(lat))
}

The above steps were taken in order to prepare the SST data for correlation and regression analysis with the three oscillations of interest. The following method employed for each of the indices is very similar to that explained above. As such, in the interest of simplicity, commentary will only accompany those steps which differ from the HADISST data above.

Employing the same process for the Atlantic Multi-decadal Oscillation (AMO), Pacific Decadal Oscillation (PDO), and El Niño-Southern Oscillation (ENSO) data.

AMO

For the climate indice data, the data was identified as being given in “months since” 1854. As such, this would have to be accounted for when formatting the date:

A.time<-ncvar_get(AMO, 'time')
A.tunits<-ncatt_get(AMO, "time",attname="units")
A.tustr<-strsplit(A.tunits$value, " ")
A.date <- as.character(as.Date(A.time*365.25/12, origin=unlist(A.tustr)[3]))

Extracting the SST data.

A.sst <- ncvar_get(AMO, 'AMO')

In this case, the print function also showed the missing value code to differ:

A.fillvalue <- ncatt_get(AMO,"AMO","_FillValue")  
A.sst[A.sst==A.fillvalue$value] <- NA 
A.missvalue <- ncatt_get(AMO,"AMO","missing_value")
A.sst[A.sst==A.missvalue$value] <- NA 
A.sst[A.sst==3000000000000000000000000000000000] <- NA

As the spatial bounds of the oscillation are not global, the global mean step was removed, but annual means still created:

A.year <- format(as.Date(A.date, format="%Y-%m-%d"),"%Y")
A.annmean <- aggregate(A.sst, by=list(A.year),FUN=mean,na.rm=TRUE)

By plotting the result, the extent of the data-set became obvious:

plot(A.annmean$Group.1, A.annmean$x, type = "l", xlab = "Time", ylab="Annual Mean SST anomaly", main = "AMO", col = "blue")

At this point, the annual means for the global SST data from the HADISST data and the annual means for the AMO could be merged in order to compare the two:

A.Merge <- merge(H.annmean, A.annmean, by = "Group.1")
colnames(A.Merge) <- c("Year", "HAD_AM", "AMV_AM")

Within this new merged dataframe, the time period for which both datasets had data present could be visually inspected in R Studio, and the result lead to the extraction of the overlapping time period (1879 - 2017):

A.Merge <- A.Merge[10:148,]

Then, the HADISST data could also be subset in order to provide an identical timeframe for which the correlation and regression results could be calculated and mapped:

A.asst <- H.asst[,,10:148]

The timeframe of interest was then

plot(A.Merge$Year, A.Merge$AMV_AM, type='l',xlab='Year',ylab='SST Anomaly',main='AMO 1879-2017', col = "blue")

At this point, the correlation analysis could be conducted. Below, the code used to create the matrices from which correlation could be compared and assessed for significance at the 0.01 alpha threshold are documented. The grid used to map the results was given the correlation and p-values as additional columns within the dataframe.

c.matrix <- matrix(NA,length(lon),length(lat))
t.matrix <- matrix(NA,length(lon),length(lat))

for (i in 1:dim(lon)) {
  for (j in 1:dim(lat)) {
    if (length(A.asst[i,j,][!is.na(A.asst[i,j,])])>2){
      c.matrix[i,j] <- cor(A.asst[i,j,], A.Merge$AMV_AM)
      p.vals <- cor.test(A.asst[i,j,], A.Merge$AMV_AM)
      t.matrix[i, j] <- p.vals$p.value
    }
  }
}

grid$corr <- as.vector(c.matrix)
grid$pval <- as.vector(t.matrix)

sig <- subset(grid[, c(1, 2, 5, 6)], pval < 0.01) 
sig <- SpatialPointsDataFrame(coords = sig[, c(1, 2)], data = sig)

The levelplot function was then used to visualise the results.

AMO.SSTCorrMap <- levelplot(corr~x*y, data=grid, xlim=c(-180,180),ylim=c(-90,90),
          col.regions = pal(100),xlab='Longitude',ylab='Latitude', main = "Correlation of AMO and HAD Annual Average SSTs") + 
  layer(sp.lines(world.coast)) +
  layer(sp.points(sig, pch = 1, cex = 0.0001, col = "black"))

AMO.SSTCorrMap

‘AMO.CorrMap’ was created. (Note: the resulting patterns will be discussed in the results section below.)

Lastly, similar steps were taken to assess significant results in the regression analysis patterns:

r.matrix <- matrix(NA,length(lon),length(lat))
s.matrix <- matrix(NA,length(lon),length(lat))

for (i in 1:dim(lon)) {
  for (j in 1:dim(lat)) {
    if (length(A.asst[i,j,][!is.na(A.asst[i,j,])])>2){
      r.lm <- lm(A.asst[i,j,]~A.Merge$AMV_AM)
      r.matrix[i,j] <- r.lm$coefficients[2]
      smm<-summary(r.lm)
      s.matrix[i, j] <- smm$coefficients[8]
    }
  }
}

Adding the results to the same data frame.

grid$reg <- as.vector(r.matrix)
grid$sig <- as.vector(s.matrix)

Creating the significance bounds.

sig <- subset(grid[, c(1, 2, 5, 6)], pval < 0.01) 
sig <- SpatialPointsDataFrame(coords = sig[, c(1, 2)], data = sig)

Plotting the result:

AMO.SSTRegMap <- levelplot(reg~x*y, data=grid, at=c(-30:30)/10,
          col.regions = pal(100),xlab='Longitude',ylab='Latitude',
          main="Regression of AMO and HAD Annual Average SSTs") + 
  layer(sp.lines(world.coast)) +
  layer(sp.points(sig, pch = 20, cex = 0.0001, col = "black"))

AMO.SSTRegMap

The same process was undertaken for PDO and ENSO:

PDO

PD.time<-ncvar_get(PDO, 'time')
PD.tunits<-ncatt_get(PDO, "time",attname="units")
PD.tustr<-strsplit(PD.tunits$value, " ")
PD.date <- as.character(as.Date(PD.time*365.25/12,origin=unlist(PD.tustr)[3]))

PD.sst <- ncvar_get(PDO, 'index')

PD.fillvalue <- ncatt_get(PDO,"index","_FillValue")  
PD.sst[PD.sst==PD.fillvalue$value] <- NA 
PD.missvalue <- ncatt_get(PDO,"index","missing_value")
PD.sst[PD.sst==PD.missvalue$value] <- NA 
PD.sst[PD.sst==-3000000000000000000000000000000000] <- NA

PD.year <- format(as.Date(PD.date, format="%Y-%m-%d"),"%Y")
PD.annmean <- aggregate(PD.sst, by=list(PD.year),FUN=mean,na.rm=TRUE)

plot(PD.annmean$Group.1, PD.annmean$x, type = "l", xlab = "Year", ylab = "SST Anomaly", main = "PDO", col = "orange")

PD.Merge <- merge(H.annmean, PD.annmean, by = "Group.1")
colnames(PD.Merge) <- c("Year", "HAD_AM", "PD_AM")

PD.Merge <- PD.Merge[16:134,]
PD.asst <- H.asst[,,30:148]

plot(PD.Merge$Year, PD.Merge$PD_AM, type='l',xlab='Year',ylab='SST Anomaly',main='PDO', col = "orange")

Calculating the correlation.

c.matrix <- matrix(NA,length(lon),length(lat))
t.matrix <- matrix(NA,length(lon),length(lat))

for (i in 1:dim(lon)) {
  for (j in 1:dim(lat)) {
    if (length(PD.asst[i,j,][!is.na(PD.asst[i,j,])])>2){
      c.matrix[i,j] <- cor(PD.asst[i,j,], PD.Merge$PD_AM)
      p.vals <- cor.test(PD.asst[i,j,], PD.Merge$PD_AM)
      t.matrix[i, j] <- p.vals$p.value
    }
  }
}

grid$corr <- as.vector(c.matrix)
grid$pval <- as.vector(t.matrix)

sig <- subset(grid[, c(1, 2, 5, 6)], pval < 0.01) 
sig <- SpatialPointsDataFrame(coords = sig[, c(1, 2)], data = sig)

PDO.SSTCorrMap <- levelplot(corr~x*y, data=grid, xlim=c(-180,180),ylim=c(-90,90),
          col.regions = pal(100),xlab='Longitude',ylab='Latitude',main="Correlation of PDO and HAD Annual Average SSTs") + 
  layer(sp.lines(world.coast)) +
  layer(sp.points(sig, pch = 1, cex = 0.0001, col = "black"))

PDO.SSTCorrMap

Calculating the regression.

r.matrix <- matrix(NA,length(lon),length(lat))
s.matrix <- matrix(NA,length(lon),length(lat))

for (i in 1:dim(lon)) {
  for (j in 1:dim(lat)) {
    if (length(PD.asst[i,j,][!is.na(PD.asst[i,j,])])>2){
      r.lm <- lm(PD.asst[i,j,]~PD.Merge$PD_AM)
      r.matrix[i,j] <- r.lm$coefficients[2]
      smm<-summary(r.lm)
      s.matrix[i, j] <- smm$coefficients[8]
    }
  }
}

grid$reg <- as.vector(r.matrix)
grid$sig <- as.vector(s.matrix)

sig <- subset(grid[, c(1, 2, 5, 6)], pval < 0.01) 
sig <- SpatialPointsDataFrame(coords = sig[, c(1, 2)], data = sig)

PDO.SSTRegMap <- levelplot(reg~x*y, data=grid, at=c(-30:30)/10,
          col.regions = pal(100),xlab='Longitude',ylab='Latitude',
          main="Regression of PDO and HAD Annual Average SSTs") + 
  layer(sp.lines(world.coast)) +
  layer(sp.points(sig, pch = 20, cex = 0.0001, col = "black"))

PDO.SSTRegMap

ENSO

EN.time<-ncvar_get(ENSO, 'time')
EN.tunits<-ncatt_get(ENSO, "time",attname="units")
EN.tustr<-strsplit(EN.tunits$value, " ")
EN.date <- as.character(as.Date(EN.time*365.25/12,origin=unlist(EN.tustr)[3]))

EN.sst <- ncvar_get(ENSO, 'Nino4')

EN.fillvalue <- ncatt_get(ENSO,"Nino4","_FillValue") 
EN.sst[EN.sst==EN.fillvalue$value] <- NA 
EN.missvalue <- ncatt_get(ENSO,"Nino4","missing_value")
EN.sst[EN.sst==EN.missvalue$value] <- NA 
EN.sst[EN.sst==-3000000000000000000000000000000000] <- NA

EN.year <- format(as.Date(EN.date, format="%Y-%m-%d"),"%Y")
EN.annmean <- aggregate(EN.sst, by=list(EN.year),FUN=mean,na.rm=TRUE)

plot(EN.annmean$Group.1, EN.annmean$x, type = "l", main = "ENSO", col = "red")

EN.Merge <- merge(H.annmean, EN.annmean, by = "Group.1")
colnames(EN.Merge) <- c("Year", "HAD_AM", "ENSO_AM")

plot(EN.Merge$Year, EN.Merge$ENSO_AM, type='l',xlab='Year',ylab='SST Anomaly',main='ENSO', col = "red")

EN.asst <- H.asst

Assessing correlation:

c.matrix <- matrix(NA,length(lon),length(lat))
t.matrix <- matrix(NA,length(lon),length(lat))

for (i in 1:dim(lon)) {
  for (j in 1:dim(lat)) {
    if (length(EN.asst[i,j,][!is.na(EN.asst[i,j,])])>2){
      c.matrix[i,j] <- cor(EN.asst[i,j,], EN.Merge$ENSO_AM)
      p.vals <- cor.test(EN.asst[i,j,], EN.Merge$ENSO_AM)
      t.matrix[i, j] <- p.vals$p.value
    }
  }
}

grid$corr <- as.vector(c.matrix)
grid$pval <- as.vector(t.matrix)

sig <- subset(grid[, c(1, 2, 5, 6)], pval < 0.01) 
sig <- SpatialPointsDataFrame(coords = sig[, c(1, 2)], data = sig)

En.SSTCorrMap <- levelplot(corr~x*y, data=grid, xlim=c(-180,180),ylim=c(-90,90),  # at=c(-1:1),
          col.regions = pal(100),xlab='Longitude',ylab='Latitude',main="Correlation of ENSO and SSTs") + 
  layer(sp.lines(world.coast)) +
  layer(sp.points(sig, pch = 1, cex = 0.0001, col = "black"))

En.SSTCorrMap

Calculating regression:

r.matrix <- matrix(NA,length(lon),length(lat))
s.matrix <- matrix(NA,length(lon),length(lat))

for (i in 1:dim(lon)) {
  for (j in 1:dim(lat)) {
    if (length(EN.asst[i,j,][!is.na(EN.asst[i,j,])])>2){
      r.lm <- lm(EN.asst[i,j,]~EN.Merge$ENSO_AM)
      r.matrix[i,j] <- r.lm$coefficients[2]
      smm<-summary(r.lm)
      s.matrix[i, j] <- smm$coefficients[8]
    }
  }
}

grid$reg <- as.vector(r.matrix)
grid$sig <- as.vector(s.matrix)

sig <- subset(grid[, c(1, 2, 5, 6)], pval < 0.01) 
sig <- SpatialPointsDataFrame(coords = sig[, c(1, 2)], data = sig)

En.SSTRegMap <- levelplot(reg~x*y, data=grid,
          col.regions = pal(100),xlab='Longitude', ylab='Latitude',
          main="Regression of ENSO and SSTs") + 
  layer(sp.lines(world.coast)) +
  layer(sp.points(sig, pch = 20, cex = 0.0001, col = "black"))

En.SSTRegMap

Sea Level Pressure (SLP)

In the interest of avoiding repetition, the SLP data was also accessed from an nc file in the same manner as the SST data above. The same process was conducted in order to extract the annually averaged SLP below:

print(SLP)

## File F:/Assignment 2 rev01/Data/slp.mon.mean.nc (NC_FORMAT_NETCDF4_CLASSIC):
## 
##      1 variables (excluding dimension variables):
##         float slp[lon,lat,time]   
##             long_name: Sea Level Pressure
##             valid_range: 870
##              valid_range: 1150
##             units: millibars
##             add_offset: 0
##             scale_factor: 1
##             missing_value: -9.96920996838687e+36
##             precision: 1
##             least_significant_digit: 1
##             var_desc: Sea Level Pressure
##             level_desc: Sea Level
##             statistic: Mean
##             parent_stat: Other
##             dataset: NCEP Reanalysis Derived Products
##             actual_range: 955.560852050781
##              actual_range: 1082.55822753906
## 
##      3 dimensions:
##         lat  Size:73
##             units: degrees_north
##             actual_range: 90
##              actual_range: -90
##             long_name: Latitude
##             standard_name: latitude
##             axis: Y
##         lon  Size:144
##             units: degrees_east
##             long_name: Longitude
##             actual_range: 0
##              actual_range: 357.5
##             standard_name: longitude
##             axis: X
##         time  Size:844   *** is unlimited ***
##             long_name: Time
##             delta_t: 0000-01-00 00:00:00
##             prev_avg_period: 0000-00-01 00:00:00
##             standard_name: time
##             axis: T
##             units: hours since 1800-01-01 00:00:0.0
##             actual_range: 1297320
##              actual_range: 1913112
## 
##     8 global attributes:
##         description: Data is from NMC initialized reanalysis
## (4x/day).  These are the 0.9950 sigma level values.
##         platform: Model
##         Conventions: COARDS
##         NCO: 20121012
##         history: Thu May  4 18:12:35 2000: ncrcat -d time,0,622 /Datasets/ncep.reanalysis.derived/surface/slp.mon.mean.nc ./surface/slp.mon.mean.nc
## Mon Jul  5 23:22:35 1999: ncrcat slp.mon.mean.nc /Datasets/ncep.reanalysis.derived/surface/slp.mon.mean.nc /dm/dmwork/nmc.rean.ingest/combinedMMs/slp.mon.mean.nc
## /home/hoop/crdc/cpreanjuke2farm/cpreanjuke2farm Thu Oct 26 23:42:16 1995 from pre.sig995.85.nc
## created 95/02/06 by Hoop (netCDF2.3)
## Converted to chunked, deflated non-packed NetCDF4 2014/09
##         title: monthly mean slp from the NCEP Reanalysis
##         References: http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanalysis.derived.html
##         dataset_title: NCEP-NCAR Reanalysis 1

SLP[["dim"]][["lat"]] 

## $name
## [1] "lat"
## 
## $len
## [1] 73
## 
## $unlim
## [1] FALSE
## 
## $group_index
## [1] 1
## 
## $group_id
## [1] 327680
## 
## $id
## [1] 0
## 
## $dimvarid
## $id
## [1] 0
## 
## $group_index
## [1] 1
## 
## $group_id
## [1] 327680
## 
## $list_index
## [1] -1
## 
## $isdimvar
## [1] TRUE
## 
## attr(,"class")
## [1] "ncid4"
## 
## $units
## [1] "degrees_north"
## 
## $vals
##  [1]  90.0  87.5  85.0  82.5  80.0  77.5  75.0  72.5  70.0  67.5  65.0
## [12]  62.5  60.0  57.5  55.0  52.5  50.0  47.5  45.0  42.5  40.0  37.5
## [23]  35.0  32.5  30.0  27.5  25.0  22.5  20.0  17.5  15.0  12.5  10.0
## [34]   7.5   5.0   2.5   0.0  -2.5  -5.0  -7.5 -10.0 -12.5 -15.0 -17.5
## [45] -20.0 -22.5 -25.0 -27.5 -30.0 -32.5 -35.0 -37.5 -40.0 -42.5 -45.0
## [56] -47.5 -50.0 -52.5 -55.0 -57.5 -60.0 -62.5 -65.0 -67.5 -70.0 -72.5
## [67] -75.0 -77.5 -80.0 -82.5 -85.0 -87.5 -90.0
## 
## $create_dimvar
## [1] TRUE
## 
## attr(,"class")
## [1] "ncdim4"

SLP[["dim"]][["lon"]][["vals"]]  

##   [1]   0.0   2.5   5.0   7.5  10.0  12.5  15.0  17.5  20.0  22.5  25.0
##  [12]  27.5  30.0  32.5  35.0  37.5  40.0  42.5  45.0  47.5  50.0  52.5
##  [23]  55.0  57.5  60.0  62.5  65.0  67.5  70.0  72.5  75.0  77.5  80.0
##  [34]  82.5  85.0  87.5  90.0  92.5  95.0  97.5 100.0 102.5 105.0 107.5
##  [45] 110.0 112.5 115.0 117.5 120.0 122.5 125.0 127.5 130.0 132.5 135.0
##  [56] 137.5 140.0 142.5 145.0 147.5 150.0 152.5 155.0 157.5 160.0 162.5
##  [67] 165.0 167.5 170.0 172.5 175.0 177.5 180.0 182.5 185.0 187.5 190.0
##  [78] 192.5 195.0 197.5 200.0 202.5 205.0 207.5 210.0 212.5 215.0 217.5
##  [89] 220.0 222.5 225.0 227.5 230.0 232.5 235.0 237.5 240.0 242.5 245.0
## [100] 247.5 250.0 252.5 255.0 257.5 260.0 262.5 265.0 267.5 270.0 272.5
## [111] 275.0 277.5 280.0 282.5 285.0 287.5 290.0 292.5 295.0 297.5 300.0
## [122] 302.5 305.0 307.5 310.0 312.5 315.0 317.5 320.0 322.5 325.0 327.5
## [133] 330.0 332.5 335.0 337.5 340.0 342.5 345.0 347.5 350.0 352.5 355.0
## [144] 357.5

lat <- ncvar_get(SLP,"lat") 
lon  <- ncvar_get(SLP,"lon")
lon <- ifelse(lon>180, lon-360,lon)
lon <- rev(lon)

Above, the SLP nc file was inspected to see if there were any differences between the formatting of it and SST. Most notably, the time units were given in hours since 1800, and the longitude was given on a 0 to 360 scale. Both of these factors would need to be standardised to match the world.coast shapefile projection, as well as the calculation of annually averaged SLP.

m.slp <- ncvar_get(SLP, "slp")

slp.date <- ncvar_get(SLP, "time")
slp.tunits<-ncatt_get(SLP,"time",attname="units")
slp.tustr<-strsplit(slp.tunits$value, " ")
slp.date<-as.character(as.Date(slp.date/24,origin=unlist(slp.tustr)[3]))

The missing value fill value was identified:

fillvalue <- ncatt_get(SLP,"slp","_FillValue")
m.slp[m.slp==fillvalue$value] <- NA
missvalue <- ncatt_get(SLP,"slp","missing_value")
m.slp[m.slp==missvalue$value] <- NA
m.slp[m.slp==-9.96921e+36] <- NA

slp.year <- format(as.Date(slp.date, format="%Y-%m-%d"),"%Y")
slp.gmean <- colMeans(m.slp, na.rm = TRUE, dims=2)
slp.annmean <- aggregate(slp.gmean,by=list(slp.year),FUN=mean,na.rm=TRUE)

slp.avslp = rowMeans(m.slp,dims=2)

levelplot(slp.avslp,col.regions = pal(100));

grid <- expand.grid(x=lon, y=lat)
grid$slp.avslp <- as.vector(slp.avslp)

levelplot(slp.avslp~x*y,grid,col.regions = pal(100),
          xlab='Longitude',ylab='Latitude',main='Average SLP') + 
  layer(sp.lines(world.coast))

slp.yrs <- slp.annmean$Group.1 
slp.nyr <- length(slp.yrs)

slp.asst <- array(NA,c(dim(lon),dim(lat),slp.nyr)) 
for (k in 1:slp.nyr) {
  slp.asst[,,k] <- rowMeans(m.slp[,,slp.year==slp.yrs[k]],na.rm=FALSE,dims=2)
}

grid$slp.an_avslp <- as.vector(rowMeans(slp.asst,na.rm=FALSE,dims=2))

Av.SLPMap <- levelplot(slp.an_avslp~x*y, data=grid,col.regions = pal(100),
          xlab='Longitude',ylab='Latitude',main='Annually Averaged SLP') + 
  layer(sp.lines(world.coast)) 

Av.SLPMap

slp.gmean <- colMeans(slp.asst, na.rm = TRUE, dims=2)
for (k in 1:slp.nyr){
  slp.asst[,,k]<-slp.asst[,,k]-matrix(slp.gmean[k],length(lon),length(lat))
}

Employing the same process for the Atlantic Multi-decadal Oscillation (AMO), Pacific Decadal Oscillation (PDO), and El Niño-Southern Oscillation (ENSO) data.

AMO and SLP

Slp.AMerge <- merge(slp.annmean, A.annmean, by = "Group.1")

As the timespan of both data-set matched, there was no need to subset either:

slp.asst.1 <- slp.asst

plot(Slp.AMerge$Group.1, Slp.AMerge$x.y, type = 'l', xlab = 'Year', ylab = 'SLP Anomaly', main = 'AMO', col = "darkgreen")

c.matrix <- matrix(NA,length(lon),length(lat))
t.matrix <- matrix(NA,length(lon),length(lat))

for (i in 1:dim(lon)) {
  for (j in 1:dim(lat)) {
    if (length(slp.asst.1[i,j,][!is.na(slp.asst.1[i,j,])])>2){
      c.matrix[i,j] <- cor(slp.asst.1[i,j,], Slp.AMerge$x.y)
      p.vals <- cor.test(slp.asst.1[i,j,], Slp.AMerge$x.y)
      t.matrix[i, j] <- p.vals$p.value
    }
  }
}

grid$corr <- as.vector(c.matrix)
grid$pval <- as.vector(t.matrix)

sig <- subset(grid[, c(1, 2, 5, 6)], pval < 0.01) 
sig <- SpatialPointsDataFrame(coords = sig[, c(1, 2)], data = sig)

AMO.SLPCorrMap <- levelplot(corr~x*y, data=grid , xlim=c(-180,180),ylim=c(-90,90),  # at=c(-1:1),
          col.regions = pal(100),xlab='Longitude',ylab='Latitude',main="Correlation of AMO and SLP") + 
  layer(sp.lines(world.coast)) +
  layer(sp.points(sig, pch = 20, cex = 0.0001, col = "black"))

AMO.SLPCorrMap

For regression:

r.matrix <- matrix(NA,length(lon),length(lat))
s.matrix <- matrix(NA,length(lon),length(lat))


for (i in 1:dim(lon)) {
  for (j in 1:dim(lat)) {
    if (length(slp.asst.1[i,j,][!is.na(slp.asst.1[i,j,])])>2){
      r.lm <- lm(slp.asst.1[i,j,]~Slp.AMerge$x.y)
      r.matrix[i,j] <- r.lm$coefficients[2]
      smm<-summary(r.lm)
      s.matrix[i, j] <- smm$coefficients[8]
    }
  }
}

grid$reg <- as.vector(r.matrix)
grid$sig <- as.vector(s.matrix)

sig <- subset(grid[, c(1, 2, 5, 6)], pval < 0.01) # 1 = lat, 2 = lon, 5 = corr, 6 = pval
sig <- SpatialPointsDataFrame(coords = sig[, c(1, 2)], data = sig)

AMO.SLPRegMap <- levelplot(reg~x*y, data=grid,
          col.regions = pal(100),xlab='Longitude',ylab='Latitude',
          main='Regression of AMO and SLP') + 
  layer(sp.lines(world.coast)) +
  layer(sp.points(sig, pch = 20, cex = 0.001, col = "black"))

AMO.SLPRegMap

PDO and SLP

Slp.PDMerge <- merge(slp.annmean, PD.annmean, by = "Group.1")
slp.asst.2 <- slp.asst

plot(Slp.PDMerge$Group.1, Slp.PDMerge$x.y, type = 'l', xlab = 'Year', ylab = 'SLP Anomaly', main = 'PDO')

c.matrix <- matrix(NA,length(lon),length(lat))
t.matrix <- matrix(NA,length(lon),length(lat))

for (i in 1:dim(lon)) {
  for (j in 1:dim(lat)) {
    if (length(slp.asst.2[i,j,][!is.na(slp.asst.2[i,j,])])>2){
      c.matrix[i,j] <- cor(slp.asst.2[i,j,], Slp.PDMerge$x.y)
      p.vals <- cor.test(slp.asst.2[i,j,], Slp.PDMerge$x.y)
      t.matrix[i, j] <- p.vals$p.value
    }
  }
}

grid$corr <- as.vector(c.matrix)
grid$pval <- as.vector(t.matrix)

sig <- subset(grid[, c(1, 2, 5, 6)], pval < 0.01) 
sig <- SpatialPointsDataFrame(coords = sig[, c(1, 2)], data = sig)

PDO.CorrMap <- levelplot(corr~x*y, data=grid , xlim=c(-180,180),ylim=c(-90,90),  
          col.regions = pal(100),xlab='Longitude',ylab='Latitude',main="Correlation of PDO and SLP") + 
  layer(sp.lines(world.coast)) +
  layer(sp.points(sig, pch = 20, cex = 0.0001, col = "black"))

PDO.CorrMap

For regression:

r.matrix <- matrix(NA,length(lon),length(lat))
s.matrix <- matrix(NA,length(lon),length(lat))


for (i in 1:dim(lon)) {
  for (j in 1:dim(lat)) {
    if (length(slp.asst.2[i,j,][!is.na(slp.asst.2[i,j,])])>2){
      r.lm <- lm(slp.asst.2[i,j,]~Slp.PDMerge$x.y)
      r.matrix[i,j] <- r.lm$coefficients[2]
      smm<-summary(r.lm)
      s.matrix[i, j] <- smm$coefficients[8]
    }
  }
}

grid$reg <- as.vector(r.matrix)
grid$sig <- as.vector(s.matrix)
sig <- subset(grid[, c(1, 2, 5, 6)], pval < 0.01)
sig <- SpatialPointsDataFrame(coords = sig[, c(1, 2)], data = sig)

PDO.RegMap <- levelplot(reg~x*y, data=grid, at=c(-30:30)/10,
          col.regions = pal(100),xlab='Longitude',ylab='Latitude',
          main='Regression of SLP and AMO') + 
  layer(sp.lines(world.coast)) +
  layer(sp.points(sig, pch = 20, cex = 0.001, col = "black"))

PDO.RegMap

ENSO and SLP

Slp.EnMerge <- merge(slp.annmean, EN.annmean, by = "Group.1")

slp.asst.3 <- slp.asst

plot(Slp.EnMerge$Group.1, Slp.EnMerge$x.y, type = 'l', xlab = 'Year', ylab = 'SLP Anomaly', main = 'ENSO')

c.matrix <- matrix(NA,length(lon),length(lat))
t.matrix <- matrix(NA,length(lon),length(lat))

for (i in 1:dim(lon)) {
  for (j in 1:dim(lat)) {
    if (length(slp.asst.3[i,j,][!is.na(slp.asst.3[i,j,])])>2){
      c.matrix[i,j] <- cor(slp.asst.3[i,j,], Slp.EnMerge$x.y)
      p.vals <- cor.test(slp.asst.3[i,j,], Slp.EnMerge$x.y)
      t.matrix[i, j] <- p.vals$p.value
    }
  }
}

grid$corr <- as.vector(c.matrix)
grid$pval <- as.vector(t.matrix)

sig <- subset(grid[, c(1, 2, 5, 6)], pval < 0.01) 
sig <- SpatialPointsDataFrame(coords = sig[, c(1, 2)], data = sig)

En.CorrMap <- levelplot(corr~x*y, data=grid , xlim=c(-180,180),ylim=c(-90,90),  # at=c(-1:1),
                         col.regions = pal(100),xlab='Longitude',ylab='Latitude',main="Correlation of ENSO and SLP") + 
  layer(sp.lines(world.coast)) +
  layer(sp.points(sig, pch = 20, cex = 0.0001, col = "black"))

En.CorrMap

r.matrix <- matrix(NA,length(lon),length(lat))
s.matrix <- matrix(NA,length(lon),length(lat))


for (i in 1:dim(lon)) {
  for (j in 1:dim(lat)) {
    if (length(slp.asst.3[i,j,][!is.na(slp.asst.3[i,j,])])>2){
      r.lm <- lm(slp.asst.3[i,j,]~Slp.EnMerge$x.y)
      r.matrix[i,j] <- r.lm$coefficients[2]
      smm<-summary(r.lm)
      s.matrix[i, j] <- smm$coefficients[8]
    }
  }
}

grid$reg <- as.vector(r.matrix)
grid$sig <- as.vector(s.matrix)
sig <- subset(grid[, c(1, 2, 5, 6)], pval < 0.01) 
sig <- SpatialPointsDataFrame(coords = sig[, c(1, 2)], data = sig)

EN.SLPRegMap <- levelplot(reg~x*y, data=grid, at=c(-30:30)/10,
          col.regions = pal(100),xlab='Longitude',ylab='Latitude',
          main='Regression of ENSO and SLP') + 
  layer(sp.lines(world.coast)) +
  layer(sp.points(sig, pch = 20, cex = 0.001, col = "black"))

EN.SLPRegMap

Assignment question: ’Indicate the significance of these patterns and include a discussion of these patterns and the role that ocean circulation may have played in producing these patterns with reference to literature.

AMO

While the contribution of ocean circulation to the dynamics of AMO processes has espoused recent debate, such as that of the opposing propositions of Clement et al. (2015) and Zhang et al. (2015), the correlation and regression patterns which emerge in these results can be thought confidently as being the result of circulation responses to the forcing of NAO manifesting in the form of the AMO (McCarthy et al. 2015).

In the AMO-SLP correlation map, the most prominent negative correlations appear over the Sahel and Middle East. In previous studies, the influence of this variability on rainfall in those regions has been identified. O’Reilly et al. (2017) refer to the observed increase in the rainy season precipitation during warm AMO phases, whereas decreases in rainfall are often associated with higher pressures over the tropical north Atlantic (Hastenrath, 1990). As this is a region with developing countries reliant on agricultural livelihoods, the significance of these atmosphere-ocean interactions manifest in socio-ecological systems in ways which can act to increase vulnerability.

For correlation with SST, a dominant pattern emerges in the Atlantic. Perhaps this shows the integration of the AMOC and the response of ocean heat transport to strong phases on the NAO. For example, such circumstances may manifest in much colder winters following a previous negative NAO phase due to the reemergence of anomalies later in the AMO (McCarthy, 2019).

There is also distinct significant correlation and regression patterns in the mid-Atlantic and Mediterranean Seas when it comes to SSTs. Delworth et al. (2016) show what could be considered as cascading impacts of the NAO cycles on AMOC and the role this may have in contributing to the longer term anthropogenic forced trend, as well as processes such as cyclone formation. In fact, Zhang et al. (2017) show the AMM (also associated with North Atlantic SSTs) being influential in the spatial distribution of North Pacific tropical cyclone-genesis; illustrating the large patterns of significant and correlation for this oscillation.

ENSO

The El-Nino Southern Oscillation is perhaps the most familiar of the indices to many due to its frequent featuring in the media, especially when it comes to reporting record breaking hot years or extremes, such as 2015 at the time.

With the results from ENSO, the correlation patterns which emerge for SSTs and SLP globally illustrate the far reaching teleconnections it can exhibit in the climate system (as is visible across the whole Pacific in the maps). There is a clear positive correlation with significance at the 0.01 level across the majority of the Pacific Ocean. This is not surprising given the nature of El-Nino and La-Nina associated events, by which trade winds and water up-welling act to affect ocean heat distribution throughout the Pacific.

The regression patterns which emerge are of great interest, especially in the context of results such as Chand et al. (2016) and their study on the effects of El-Nino in the future. In their multi-modelling experiment, the results show an increase in tropical cyclone frequency during El-Nino events compared to now around many of the small Pacific Island states such as Hawaii. This, given the strong positive correlation obvious with SLP and SSTs, shows the significance of knowing these patterns and the potential benefit this information can have for adaptation planning in these small vulnerable island nations.

Lastly, Jiménez-Muñoz et al. (2016) bring another important perspective to the forefront regarding the pervasive facts in the Southern Oscillation: Amazonian drought. With the correlation and regression maps, the pattern of significance over the northern part of South America is evident and perhaps could be one of the processes contributing to these conditions over time. With the issue of Amazon Forest Dieback being prominent in research around tipping points in the climate system and abrupt change (Mahli et al., 2008; Zemp et al., 2017), the integration of these ocean systems with terrestrial ecosystems even becomes evident, showing that the spatial distribution of impacts is hugely significant.

PDO

As Newman et al. (2015) point out, the PDO is a result of a multitude of processes, with the Kurishio-Oyashio associated rossby wave propagation, as well as fluxes and Ekman transport related to the Aleutian low being two main contributors resulting associated ocean circulation.

In the correlation maps for SST and SLP, the sheer scale of the significant positive zonal influence on SLP is evident at the equator and the tropics. Interestingly, a strong negative correlation is present for SLP off Japan. When it comes to SSTs, correlation and regression patterns of significance are evident throughout the Pacific as well as the Indian Ocean - perhaps linked to Meridional overturning circulation globally.

The consequence and significance of this distribution is evident in the results of Hu et al. (2017) who use Earth System Model runs to assess how these pressure driven (wind) oscillations can contribute similarly to ocean heat redistribution, but particularly in the context of the connecting of the PDO and AMO by the AMOC - perhaps explaining the presence of significant negative correlations in the Atlantic for the PDO off the coast of Ireland.

The significance of this may manifest in its impact on SST anomalies in the Pacific, Indian, or Atlantic Ocean. As Hu et al. (2017) also reference, the potential for modulation of external forcing shows the impact of this variability, with Zhou et al. (2017) attributing the long-term sub-surface cooling observed in the Indian Ocean to PDO mechanics in spite of surface temperate increases throughout the same time period.

Conclusion

In conclusion, this assignment set out to map correlation and regression patterns for three climate indices with respect to long-term SLP and SST data. The methods and code used to create the resulting outputs on R is included in this R Markdown document in the ethos of reproducibility, and to illustrate steps taken. The results of this assignment highlight the interconnected nature of ocean heat distribution and the impact of natural climate variability across all parts of the world. With correlation and regression analysis such as above, it is hoped this information can allow for better climate projections through improved understanding of the ocean-atmospheric climate system, and to a degree, the potential for predictability on a short timescale. Perhaps useful for potential climate services or to assist in adaptation planning.

Bibliography

Bader, J., & Latif, M. (2003) The impact of decadal scale Indian Ocean sea surface temperature anomalies on Sahelian rainfall and the North Atlantic Oscillation. Geophysical Research Letters [online]. 30(22). Available at: https://agupubs.onlinelibrary.wiley.com/doi/full/10.1029/2003GL018426 (accessed 24 April 2019).

Chand, S. S., Tory, K. J., Ye, H., & Walsh, K. J. (2017) Projected increase in El Nino-driven tropical cyclone frequency in the Pacific. Nature Climate Change [online]. 7(2), 123. Available at: https://nab.vu/sites/default/files/documents/Chand%20et%20al.%202017.pdf (accessed 28 April 2019).

Clement, A., Bellomo, K., Murphy, L. N., Cane, M. A., Mauritsen, T., Rädel, G., & Stevens, B. (2015) The Atlantic Multidecadal Oscillation without a role for ocean circulation. Science [online]. 350(6258), 320-324. Available at: https://science.sciencemag.org/content/sci/350/6258/320.full.pdf (accessed 26 April 2019).

Delworth, T. L., Zeng, F., Vecchi, G. A., Yang, X., Zhang, L., & Zhang, R. (2016) The North Atlantic Oscillation as a driver of rapid climate change in the Northern Hemisphere. Nature Geoscience [online]. 9(7), 509. Available at: https://www.researchgate.net/profile/Liping_Zhang19/publication/304336432_The_North_Atlantic_Oscillation_as_a_driver_of_rapid_climate_change_in_the_Northern_Hemisphere/links/5c059588a6fdcc315f9ad47c/The-North-Atlantic-Oscillation-as-a-driver-of-rapid-climate-change-in-the-Northern-Hemisphere.pdf (accessed 28 April 2019).

Hastenrath, S. (1990) Decadal-scale changes of the circulation in the tropical Atlantic sector associated with Sahel drought. International Journal of Climatology [online] 10(5), 459-472. Available at: https://rmets.onlinelibrary.wiley.com/doi/abs/10.1002/joc.3370100504 (accessed 28 April 2019).

Hu, Z., Hu, A., & Hu, Y. (2018) Contributions of Interdecadal Pacific oscillation and Atlantic multidecadal oscillation to Global Ocean heat content distribution. Journal of Climate [online]. 31(3), 1227-1244. Available at: https://journals.ametsoc.org/doi/pdf/10.1175/JCLI-D-17-0204.1 (accessed 28 April 2019).

Jiménez-Muñoz, J. C., Mattar, C., Barichivich, J., Santamaría-Artigas, A., Takahashi, K., Malhi, Y., … & Van Der Schrier, G. (2016) Record-breaking warming and extreme drought in the Amazon rainforest during the course of El Niño 2015-2016. Scientific reports [online]. 6, 33130. Available at: https://www.nature.com/articles/srep33130 (accessed 28 April 2019).

Malhi, Y., Aragão, L. E., Galbraith, D., Huntingford, C., Fisher, R., Zelazowski, P., … & Meir, P. (2009) Exploring the likelihood and mechanism of a climate-change-induced dieback of the Amazon rainforest. Proceedings of the National Academy of Sciences [online]. 106(49), 20610-20615. Available at: https://www.pnas.org/content/pnas/106/49/20610.full.pdf (accessed 29 April 2019).

McCarthy, G. D., Haigh, I. D., Hirschi, J. J. M., Grist, J. P., & Smeed, D. A. (2015) Ocean impact on decadal Atlantic climate variability revealed by sea-level observations. Nature [online]. 521(7553), 508. Available at: http://nora.nerc.ac.uk/id/eprint/510937/1/mccarthy_manuscript_post_print.pdf (accessed 26 April 2019).

Munk, W. H. (1950) On the wind-driven ocean circulation. Journal of meteorology [online]. 7(2), 80-93. Available at: https://journals.ametsoc.org/doi/pdf/10.1175/1520-0469(1950)007%3C0080:OTWDOC%3E2.0.CO%3B2 (accessed 22 April 2019).

National Oceanic and Atmospheric Administration (NOAA) (2019) Energy in the Ocean and Atmosphere [online]. Available at: https://oceanservice.noaa.gov/education/pd/oceans_weather_climate/energy_oceans_atmosphere.html (accessed 23 April 2019).

Nelson, S. A. (2016) The Ocean-Atmosphere System [online]. Available at: http://www.tulane.edu/~sanelson/Natural_Disasters/oceanatmos.htm (accessed 22 April 2019).

Newman, M., Alexander, M. A., Ault, T. R., Cobb, K. M., Deser, C., Di Lorenzo, E., … & Schneider, N. (2016) The Pacific decadal oscillation, revisited. Journal of Climate [online]. 29(12), 4399-4427. Available at: https://journals.ametsoc.org/doi/pdf/10.1175/JCLI-D-15-0508.1 (accessed 28 April 2019).

O’Reilly, C. H., Woollings, T., & Zanna, L. (2017) The dynamical influence of the Atlantic Multidecadal Oscillation on continental climate. Journal of Climate [online]. 30(18), 7213-7230. Available at: https://journals.ametsoc.org/doi/pdf/10.1175/JCLI-D-16-0345.1 (accessed 27 April 2019).

Pachauri, R. K., & Meyer, L. (2014). Climate change 2014 Synthesis Report-Summary for Policymakers [online]. Available at: https://www.ipcc.ch/site/assets/uploads/2018/02/AR5_SYR_FINAL_SPM.pdf (accessed 23 April 2019).

Peterson, R. G., & White, W. B. (1998) Slow oceanic teleconnections linking the Antarctic circumpolar wave with the tropical El Nino Southern Oscillation. Journal of Geophysical Research: Oceans [online]. 103(C11), 24573-24583. Available at: https://agupubs.onlinelibrary.wiley.com/doi/pdf/10.1029/98JC01947 (accessed 25 April 2019).

Zhang, W., Vecchi, G. A., Villarini, G., Murakami, H., Rosati, A., Yang, X., … & Zeng, F. (2017) Modulation of western North Pacific tropical cyclone activity by the Atlantic Meridional Mode. Climate dynamics [online]. 48(1-2), 631-647. Available at: https://link.springer.com/article/10.1007/s00382-016-3099-2 (accessed 29 April 2019).

Zhang, R., Sutton, R., Danabasoglu, G., Delworth, T. L., Kim, W. M., Robson, J., & Yeager, S. G. (2016) Comment on “The Atlantic Multidecadal Oscillation without a role for ocean circulation”. Science [online]. 352(6293), 1527-1527. Available at: https://science.sciencemag.org/content/sci/352/6293/1527.1.full.pdf (accessed 27 April 2019).

Zemp, D. C., Schleussner, C. F., Barbosa, H. M., Hirota, M., Montade, V., Sampaio, G., … & Rammig, A. (2017) Self-amplified Amazon forest loss due to vegetation-atmosphere feedbacks. Nature communications [online]. 8, 14681. Available at: https://www.nature.com/articles/ncomms14681 (accessed 29 April 2019).

Zhou, X., Alves, O., Marsland, S. J., Bi, D., & Hirst, A. C. (2017) Multi-decadal variations of the South Indian Ocean subsurface temperature influenced by Pacific Decadal Oscillation. Tellus A: Dynamic Meteorology and Oceanography [online]. 69(1), 1308055. Available at: https://www.tandfonline.com/doi/full/10.1080/16000870.2017.1308055 (accessed 29 April 2019).