Supplementary material for ‘The Little Ice Age climate of New Zealand reconstructed from Southern Alps cirque glaciers: a synoptic type approach” by Lorrey et al.

S1. Sources of geospatial data and precision of Google Earth imagery

Sources and precision of GE geospatial data are not provided by Google™ (including digital elevation model (DEM) details; Potere, 2008). General consensus among the GE user community is that a DEM from NASA’s Shuttle Radar Topography Mission (SRTM; Farr et al., 2007) is employed. SRTM data exist for 80% of the Earth’s surface between 60°N and 56°S. Most countries aside from the United States have a 90 m resolution. This means that for NZ there is probably 90 m horizontal spacing between elevation control points, and elevation values between those points are interpolated (Farr et al., 2007). Vertical precision in the 90 m SRTM dataset is approximately 30 m (Farr et al., 2007). The SRTM dataset has known issues related to shadowing or foreshortening caused by steep slopes that introduce voids into the dataset. For high elevation regions, vegetation snow and ice cover are known factors that complicate radar wave propagation, so it is possible that the ground surface has not always been mapped. Radar can penetrate into snow and ice, so depending on the state of the snow/ice pack, the radar may have measured elevations from the top of the pack or from the ground surface underneath, or from partway through the ice (Farr et al., 2007).

Some higher resolution DEM data has been incorporated into GE through third party involvement in Google’s Map Content Partner Program (e.g. a 40 cm gridded DEM was provided to GE by Marin County, California). It cannot be established whether the Southern Alps has the SRTM DEM (90 m) or a third-party high-resolution DEM. In addition to the ambiguity surrounding the DEMs in use in GE, there are issues surrounding misplacement of imagery in the programme. Some images have been found to be incorrectly positioned relative to adjacent images, resulting in, for example, disjointed roads and coastlines (Potere, 2008). If an image is misaligned to the underlying DEM, this creates problems for obtaining elevation data for a particular point identified on the image.

S2. Comparison of ensemble composites with and without use of weighting

The decision to use analog season repeats as a weighting mechanism for constructing the ensemble composite patterns for z1000 and SSTa was knowingly made, despite the artificial inflation of t-test significance. The argument for doing this is that 1.) use of weighting helps to indicate which anomalies were the most common between sites and 2.) it does not greatly distort the spatial patterns that result for weighted result relative to the unweighted ensemble composites. A plot of the z1000 result for the Southern Hemisphere (Figure S1), which can be compared directly with Figure 8 in the text, is used to demonstrate the similarities of the two approaches.

Figure S1. Plot of geopotential height anomaly at 1000hPa (z1000) for the Southern Hemisphere during the LIA based on analog years identified by temperature reconstructions for 22 Southern Alps glaciers. The 90th and 95th percent (dashed and solid black lines) confidence intervals are superposed on the plot. While the extent of the confidence in the plot is less than that seen in Figure 8, overall the lack of weighting for individual analog years used in this diagram results in a spatial pattern (including location of ‘highs’ and ‘lows’ and strength of those features) that is similar to what is seen in Figure 8.

The mean austral summer temperature anomaly derived from 22 Southern Alps cirque glaciers (~0.56°C +/- 0.29°C) and that from Oroko Swamp tree rings (-0.42°C; using data from Cook et al. 2002) is different from borehole reconstructions of mean annual temperature anomalies in the Eastern South Island (-0.9°C; Whiteford et al. 1996). In addition, SSTa reconstructions provided by Mann et al. (2009) suggest colder temperatures offshore of the west coast of the South Island during the LIA than the result from this study and also that of Cook et al. (2002). We have suggested that the differences between the summer reconstructions (this study; Cook et al., 2002) and the annual reconstructions (Whiteford et al. 1996; Mann et al. 2009) could be due to colder winter temperature anomalies relative to summer temperature anomalies during the LIA. We hypothesised that colder SSTs, which would have included and were perhaps been driven in-part by increased sea ice concentration and more expansive sea ice extent during winter relative to summer in regions upwind of NZ, promoted the difference between summer and annual temperature anomalies during the LIA. We have taken the opportunity to interrogate a global climate model simulation for the LIA to examine sea ice changes relative to the modern period to test this hypothesis.

The CSIRO Mk3L model indicates more extensive sea ice and increased sea ice concentrations in general to the south of Australia during the LIA (Figure S2). Winter sea ice extent is also further north and more concentration than for summer. We suggest that this seasonal contrast may have had the effect of driving down winter temperature anomalies (and therefore the overall annual average) relative to summer during the LIA, which would result in lower mean annual temperatures relative to summer as well.

Figure S2. Anomaly in sea ice concentration (percent), average for 1450-1850 AD (relative to the 1961-1990 average). (a) annual mean, (b) December-January-February (austral summer) mean, and (c) June-July-August (austral winter) mean. Values are only shown where the anomalies are greater than or equal to 15%. The values were derived from a transient simulation of the period 851 to 2000 CE conducted using version 1.2 of the CSIRO Mk3L climate system model. The experimental design follows the PMIP3/CMIP5 protocol for the Last Millennium and Historical experiments, and incorporates orbital, greenhouse gas, solar and volcanic forcings.

S.3 Caveats for calculation of palaeoequilibrium lines for the LIA without reconstructed glacier profiles

The use of cirque glacier mid-points (Chinn 1996), which is effectively a THAR method of reconstructing ELAp (Meierding, 1982), via topographic maps and Google Earth may slightly overestimate the temperature depression during the LIA. The tie point of each glacier at the cirque headwall may have been at a higher elevation during the LIA, and the midpoint along each glacier profile may also have been higher that what is shown on the maps or DEM as depicted by Google Earth. In this situation, a relative comparison between profiles with and without ice would yield ELA changes that could be considered relative (equivocal). However, a tie point positioned as it is during the present and a higher LIA profile would result in a reduced temperature depression during the LIA. However, because cirque glaciers have relatively thin ice coverage, because of the errors in the elevation models, the relatively subtle changes in ELA from present to the LIA, and the use of quintiles for defining the analog seasons rather than direct analog selection all probably serve to mask any possible deficiency related to the absence of a reconstructed ice profile (and the relative error that would be associated with reconstructed temperature based on ELAp). In addition, the comparison to independent proxy data to determine if realistic values have been reconstructed is recommended, and this step has been taken in this study to corroborate the results. That said, the application of the approach used in this study to longer time scales (i.e. outside the late Holocene) must be done cautiously. It is recommended that more elegant approaches, such as use of an accumulation area ratio (AAR) and reconstructed ice profiles, may yield more robust results for ELA-based temperature reconstructions over earth orbital scales for New Zealand.

Supplement References not in main body text

Farr, T.G., Rosen, P.A., Caro, E., Crippen, R., Duren, R., Hensley, S., Kobrick, M., Paller, M., Rodriguez, E., Roth, L., Seal, D., Shaffer, S., Shimada, J., Umland, J., Werner, M., Oskin, M., Burbank, D. and Alsdorf, D. 2007. The Shuttle Radar Topography Mission, Reviews of Geophysics, 45, RG2004, doi:10.1029/2005RG000183.

Potere, D. 2008. Horizontal positional accuracy of Google Earth’s high-resolution imagery archive, Sensors, 8: 7973-7981.