Tree-ring stable carbon isotope-based June–September maximum temperature reconstruction since AD 1788, north-west Thailand

Abstract

The first study of tree-ring stable carbon isotopes in Thailand has demonstrated that stable carbon isotope in northwestern Thailand represents a promising proxy for the temperature reconstruction of core-monsoon periods. A tree-ring δ¹³C chronology was constructed based on four cores covering the period of 1788–2013. After removing the long-term decreasing trend reflecting atmospheric CO₂ concentrations, the ∆¹³C chronology was able to capture both temperature and hydro-climate signals. ∆¹³C chronology showed particularly strong and significant negative correlation (r = –0.62, p < 0.0001) with June–September maximum temperature (CRU TS 3.24). The maximum temperature was reconstructed, which explained 37.8% of the variance in the instrumental maximum temperatures over the period of 1901–2013. The maximum temperature reconstruction revealed that four cooler and three warmer periods, as well as a slightly increasing temperature trend, occurred during the late seventeenth to mid-eighteenth centuries, which were followed by severe temperature fluctuations during the twentieth century century. While the sea surface temperature anomaly in the Indian Ocean might not affect the maximum temperature, its unstable relationship with the El Niño-Southern Oscillation (ENSO) was detected. In addition, a close relationship was observed between the maximum temperature and ENSO during the negative phase of the Pacific Decadal Oscillation (PDO), but this relationship was lost during the positive phase of the PDO.

Keywords:

• north-west Thailand
• tree-ring cellulose carbon isotope
• Pinus merkusii
• maximum temperature

1. Introduction

The spread of dense population, civilization and biodiversity throughout South-East Asia has been shaped and characterized by the influence of the Asian monsoon. Since the beginning of modern civilization, anthropogenic activity has been a potential factor in driving global climate; it tends to increase the global temperature, which causes variable changes in monsoon systems (IPCC, 2013) In Thailand, 80% of annual rainfall occurs during the summer monsoon season (May–October), when the feast or famine of a population depends on monsoonal rainfall. Contradictory changes in the Asian monsoon over Thailand under the background of global warming have been investigated due to the interaction between the climate of the Pacific and Indian Oceans and driving forces such as the El Niño-Southern Oscillation (ENSO) and the warm pool in the Indian Ocean (the Indian Ocean Dipole). Climate records obtained using instruments and natural palaeoclimate archives are the keys to understanding the spatial and temporal changes in monsoon mechanisms and their magnitudes. Studies have used instrumental climate data to reveal the relationship between ENSO and the onset of monsoons, which suggested that weak monsoons occurred in this region during El Niño events, and vice versa (Rupakumar and Pant, 1997; Hamada et al., 2002). In addition, an unsteady relationship has existed between ENSO and the Indian summer monsoon since 1980 (Kumar et al., 1999). Due to the south-eastward shift in the descending limb of the Walker circulation during El Niño events, weak correlations between ENSO and Indian summer monsoon rainfall after 1980 have also appeared in summer monsoon rainfall in northern Thailand, and its correlation with rainfall increased in the centre of Thailand (Singhrattna et al., 2005).

There are many natural archives of palaeoclimatic information, including ice cores, lake sediment (varves), tree-rings, coral reefs, and pollen, among others (Pielou, 1991; Agashe, 1995; Bradley, 1999; Fagan, 2000). Tree-ring analyses, unlike other proxy records, have two great advantages: (1) each ring yields an exactly annual resolution and (2) multiple samples enhance the signal to noise ratio. In tropical and subtropical regions, where rainfall and temperature are not limiting factors for tree growth, proxy data, such as trees that build obvious annual rings, are limited. Three tree species (Tectona grandis, Pinus merkusii and Pinus kesiya) have been used in dendrochronology studies in Thailand and have been considered to be promising climate proxies that record precipitation, temperature and Palmer Drought Severity Index (PDSI) data (Buckley, Barbetti, et al., 1995; Pumijumnong et al., 1995; Buckley, Duangsathaporn, et al. 2007; Buckley, Palakit, et al., 2007; Pumijumnong and Eckstein, 2010; Pumijumnong, 2012). However, a limitation of traditional tree-ring width analysis in this region is that these data usually only reflect climate signals during premonsoon periods.

In recent decades, using new advances in isotopic geochemistry, including a better understanding of the physiological controls on isotopic variations (McCarroll and Loader, 2004; McCarroll et al., 2009; Leavitt, 2010), an increasing number of researchers have been able to reconstruct climate factors using stable isotope ratios in tree-rings in South-East Asia. Tree-ring stable isotopes highlight two advantages: (1) the low-frequency climate signals are well preserved since detrending processes are unnecessary during their reconstructions (McCarroll and Loader, 2004; Gagen et al., 2007; Szymczak et al., 2012); (2) far less samples are required (4–5 in general) for building a reliable chronology, comparing with other tree-ring parameters (Gagen et al., 2007; Leavitt, 2010; Liu et al., 2014). Tree-ring stable oxygen isotope ratios (δ¹⁸O) have been rapidly measured throughout South-East Asia. In north-west Thailand, Pinus kasiya was used to study intra-annual variations; it exhibited negative correlations with local July–November precipitation (Zhu et al., 2012). Subsequently, July–October precipitation was reconstructed using Pinus merkusii δ¹⁸O which showed close relationship with ENSO during the period 1871–2000 (Xu et al., 2015). Sano et al. (2012) reconstructed an annual multivariate ENSO Index over the past 300 years using Fokienia hodginsii cellulose δ¹⁸O data in Vietnam, and their δ¹⁸O chronology exhibited significant correlations with precipitation and PDSI. Furthermore, Xu et al. (2011) indicated that the δ¹⁸O chronology in Vietnam showed good agreement with that in Northern Laos. The latter showed a significantly negative correlation with the May–October PDSI and revealed the significant impact of ENSO on tree-ring cellulose in Laos. Although tree ring δ¹⁸O data in South-East Asia have been proven to represent a reliable hydroclimatic archive, in order to further understand climate characteristics, stable carbon isotope ratios in tree-rings (δ¹³C) must be studied. The δ¹³C values in tree-rings are a result of the fractionation of the stable carbon isotope composition of the atmosphere due to the balance between stomata conduction, the photosynthesis rate and carboxylation (McCarroll and Loader, 2004). Therefore in many previous studies, the climatic signals that were captured by δ¹³C in tree rings include relative humidity, soil moisture, light level and; particularly, temperature (Leavitt et al., 2002; Nakatsuka et al., 2004 Gagen et al., 2007; Liu et al., 2012, 2014; Zhao et al., 2014; Bégin et al., 2015) (which are related to stomata conduction and photosynthesis).

In this article, we present a maximum temperature reconstruction spanning 228 years based on δ¹³C values in a pine tree ring in northwestern Thailand. This reconstruction captured the core-monsoon temperature and revealed increases in temperature fluctuations during the twentieth century. In addition, this 228-year-long temperature reconstruction was used to investigate the influence of the teleconnection of the Pacific and Indian Oceans on northwestern Thailand.

2. Materials and methods

2.1. Study area and climatology

Merkus pines (Pinus merkusii) were collected from the Mae Surin arboretum (18°52′34′′N, 97°56′36′′E) which is located in Mae Hong Son province (MHS), northwestern Thailand (Fig. 1). The geography of this site is characterized as a sandstone plateau (612 m a.s.l). The dominant species of the sampling site are two-needle leaf pine (Pinus merkusii) and dipterocarp (Shorea obtuse, Shorea siamensis) (Bumrungsak, 2003; Santisuk, 2013). In general, the climate of the sampling site is influenced by the south-west monsoon (summer monsoon) and the north-east monsoon (winter monsoon). The south-west monsoon appears from May to October (described as a rainy or wet monsoon season) and brings a warm, moist air mass from the Indian Ocean. The north-east monsoon appears from November to February (described as a cool season) and brings a cold dry air mass from China. According to the climate diagram obtained from the MHS meteorological station during 1951–2013 AD (Fig. 2), the dry period (in which there is less precipitation than evaporation) occurs from November to March, while the wet period (in which there is more precipitation than evaporation) occurs from May to October. The mean annual rainfall is 1228 mm; observations indicate that the highest rainfall occurs in August and the lowest rainfall occurs in February. The 63-year mean temperature is 25 °C; the highest mean temperature occurs in April (30 °C) and the lowest occurs in January (21 °C) (Thai Meteorological Department, 2015).

Fig. 1. The sampling site: Mae Surin arboretum, Mae Hong Son province, northwestern Thailand.

Fig. 2. Monthly mean temperature, precipitation and relative humidity data since 1951 in Mae Hong Son.

The samples were collected in May 2014, during the transition period between dry and wet seasons. An increment borer with a diameter of 10 mm was used to collect the samples at breast height. The cores were dried at room temperature and sanded until their cross-sections revealed boundaries of rings that were clear enough to produce tree ring measurements. The cores were measured using a LinTab measuring system (with a precision of 0.001 mm) and plotted using a TSAP-Win program (Rinn, 2011). To cross-date the tree-ring width series, each individual series was visually compared with the plots graphed using the TSAP program, which were statistically checked by the COFECHA program (Holmes, 1983; Grissino-Mayer, 2001).

2.2. Cellulose extraction and isotopic measurements

Four cores from four trees with clear ring boundaries and no missing ring were selected for isotopic analysis. During the most recent 100 years, the rings are very narrow, and they contain a few layers of late wood; although the rings during the younger period are much wider, the boundary between early wood and late wood is not clear due to the high amounts of resin in the wood. Therefore, we used the whole annual ring to measure stable carbon isotopes. Each core was cut into very small pieces (with a thickness of <0.5 mm for potential chemical reactions), year by year, using a razor blade under a stereoscope. The Jayme-Wise method was employed to extract α-cellulose (Green, 1963). The homogenized cellulose samples were wrapped in tin capsules, and their stable carbon isotope ratios were measured using a continuous flow isotope ratio mass spectrometer (Flash 2000 and Delta V Advantage). The ratios of ¹³C to ¹²C are presented using δ notation, which represents the deviation (‰) from a reference standard (IAEA CH3) based on the following equation:(1) $${\delta }^{13}C=\left(R_{\text{sample}}/R_{\text{standard}}-1\right)\times 100$$(1)

where R_sample and R_standard are the ¹³C/¹²C ratios of the sample and the standard, respectively.

After each individual δ¹³C time series was obtained, the first several decades were removed to avoid the juvenile effect (McCarroll and Loader, 2004; Leavitt, 2010). Then, the remaining individual δ¹³C time series were defined as PM1b (1788–2013 AD), PM2a (1810–2013 AD), PM4b (1788–2013 AD) and PM5a (1788–2013 AD). Because the atmospheric CO₂ concentration has continuously increased since the Industrial Revolution, ¹³C has become increasingly depleted in the atmosphere, which has caused the long-term declining trend of tree-ring δ¹³C (McCarroll and Loader, 2004; McCarroll et al., 2009). To reveal the common climatic signal in the δ¹³C time series, the effects of atmospheric CO₂ should be removed. Here, we used the discrimination method (∆¹³C, McCarroll and Loader, 2004) to calibrate each individual δ¹³C time series based on the δ¹³C atmospheric CO₂ measurements derived from ice core measurements and direct atmospheric monitoring, using the following equation:(2) $${\mathrm{\Delta }}^{13}\text{C}=\left({\delta }^{13}{\text{C}}_{\text{a}}-{\delta }^{13}{\text{C}}_{\text{tree}}\right)/\left(1-{\delta }^{13}{\text{C}}_{\text{tree}}/1000\right)$$(2)

where δ¹³C_a is the δ¹³C of the atmospheric CO₂, and δ¹³C_tree is the δ¹³C value in the tree-ring.

To diminish the differentiation of the physiological response of each tree to ∆¹³C values, each ∆¹³C series was compressed by the average value of each series to yield one single isotopic discrimination value ranging from 1788 to 2013 AD. The running expressed population signal (EPS), which is the chronology confidence estimation commonly used in dendrochronology studies and is based on the number of trees and the mean correlations between trees, was performed using a window of 30 years and a sliding step of 15 years.

2.3. Climate data and statistical analysis

To evaluate the climatic response of tree-ring δ¹³C, two sources of climatological data were used in this study. First, monthly climatological data (precipitation, temperature, relative humidity and number of rainy days) from 1951–2013 AD were obtained from the MHS meteorological station (269 m a.s.l.), which is located approximately 40 km from the sampling site. Second, six climatic parameters, namely, cloud cover, Diurnal Temperature Range (DTR), precipitation, mean temperature (T_mean), minimum temperature (T_min), and maximum temperature (T_max), were obtained from the climatic Research Unit time-series (CRU TS) 3.24 (University of East Anglia Climatic Research Unit et al., 2008) in 0.5° × 0.5° gridded data-sets (97°45′–98°15′E, 18°45′–19°15′N) from 1901–2013.

Statistical tests used to investigate the relationship between mean ∆¹³C time series and monthly climatic data were performed using the DENDROCLIM 2002 software (Biondi and Waikul, 2004). DENDROCLIM estimates significant correlation coefficient values using a bootstrapped correlation function for a single interval (i.e. an entire time series correlation) and multiple intervals (i.e. a moving interval and an evolutionary interval). In addition, monthly climatic data were combined between two different periods to determine the influence of seasonal climatic factors on mean ∆¹³C time series. Then, the climate parameter that showed the best correlation with the ∆¹³C series was reconstructed using a simple linear regression. Split-period calibration-verification was used to test the stability and reliability of the climate reconstruction model (Fritts, 1976; Fritts, 1991; Cook and Kairiukstis, 1990; Cook et al., 1999). The correlation coefficient (r), sign test, reduction of error test (or RE) (Fritts, 1976), and coefficient of efficiency (or CE) (Cook et al., 1999) were used to verify this reconstruction, along with the split-period calibration-verification.

3. Results

3.1. δ¹³C and ∆¹³C

The 4 δ¹³C series were significantly inter-correlated (p < 0.01), as they yielded r values ranging from 0.66 to 0.84 (Table 1). The mean values of the δ¹³C series fell between –23.58 and –23.48‰, and the differences between the maximum and minimum values of the time series from PM1b, PM2a, PM4b and PM5a were 3.85, 3.58, 4.23 and 3.65‰, respectively. All 4 δ¹³C series showed slightly deceasing trends throughout almost all of their series, and the linear regression of the average 4 δ¹³C series over time (1788–2013 AD) yielded a decreasing rate of 0.0046‰ yr⁻¹. However, after 1850 AD (which represents the beginning of the Industrial Revolution), the rate of decrease accelerated to 0.01‰ yr⁻¹, and it dramatically decreased to a rate of 0.028‰ yr⁻¹ after 1960. The Pearson’s correlation between the 4 δ¹³C series and the δ¹³C values of atmospheric CO₂ ranged from 0.53 to 0.76 (p < 0.01) (Fig. 3a).

Table 1. The correlation between the original series (δ¹³C) and detrended series (∆¹³C) from Mae Hong Son (1788–2013 AD).

		δ^13C				∆^13C
		PM1b	PM2a	PM4b	PM5a	PM1b	PM2a	PM4b	PM5a
δ^13C	PM1b	1
	PM2a	0.83	1
	PM4b	0.79	0.84	1
	PM5a	0.66	0.83	0.80	1
	Average	0.89	0.94	0.93	0.89
∆^13C	PM1b					1
	PM2a					0.65	1
	PM4b					0.62	0.65	1
	PM5a					0.53	0.66	0.65	1
	Average					0.86	0.86	0.85	0.83

Note: Correlation is significant at the 0.01 level.

Fig. 3. (a) Average δ¹³C (black) values from 4 trees (grey) and δ¹³C values in the atmosphere (blue); (b) ∆¹³C series after removing the long-term decreasing trend of δ¹³C in the atmosphere; (c) the EPS and Rbar statistics of ∆¹³C.

After being corrected to remove the long-term atmospheric decline in δ¹³C (Fig. 3b) using the discrimination equation (Equation (2(2) $${\mathrm{\Delta }}^{13}\text{C}=\left({\delta }^{13}{\text{C}}_{\text{a}}-{\delta }^{13}{\text{C}}_{\text{tree}}\right)/\left(1-{\delta }^{13}{\text{C}}_{\text{tree}}/1000\right)$$(2) )), the inter-correlation between each ∆¹³C series decreased to be lower than the correlation between the original δ¹³C series (i.e. from r = 0.66–0.84 and p < 0.01 to r = 0.53–0.66 and p < 0.01) (Table 1). Considering the high correlation between the original δ¹³C series and the δ¹³C of atmospheric CO₂, as well as the decrease in the inter-correlation values after removing the atmospheric decline in δ¹³C, these data indicated that the long-term declining trends in all four original δ¹³C series from the 4 Pinus merkusii trees were strongly affected by the decreasing trend of atmospheric δ¹³C due to the increasing concentrations of CO₂ in the atmosphere.

Although the inter-correlation values of the ∆¹³C series were lower than those of the original δ¹³C series, the mean inter-series correlation (Rbar) among the individual ∆¹³C series ranged from 0.52 to 0.77 over the period of 1788–2013 AD, and the running expressed population signal (EPS) values, based on a 30-year window and a sliding step of 15 years over the period of 1788–2013 AD, ranged from 0.81 to 0.93 (Fig. 3c). The analysis of Pinus merkusii and Pinus kasiya tree-ring widths from 16 sites (with approximately 26 trees per site) along the northwestern region to Upper Peninsula of Thailand showed tree ring-width correlations ranging from 0.31 to 0.50 (Pumijumnong and Eckstein, 2010). Thus, the ∆¹³C series yielded a potentially stronger correlation than the tree ring-width time series while also using fewer samples.

3.2. ∆¹³C climatic response

Four individual ∆¹³C series were combined into a master series covering a time period of 226 years (Fig. 3b). Pearson’s correlation was used to identify climatic responses (Fig. 4).

Fig. 4. Correlations between tree-ring ∆¹³C chronology and climatic factors obtained from the MHS station (mean temperature, rainy days and relative humidity) and the CRU TS3.24 (diurnal temperature and maximum temperature) at the 95% significance confidence level. The grey area is the 95% confidence limit after the Bonferroni adjustment.

During the common record period of two climatological data-sets (1951–2013 AD), the ∆¹³C chronology showed a weak positive correlation with precipitation in July (a month into the first half of monsoon season) on both a local (MHS meteorological station; r = 0.28, p < 0.05, n = 63) and regional scale (CRU TS 3.24; r = 0.28, p < 0.05, n = 63), but it exhibited higher correlations with relative humidity and rainy days (MHS meteorological station; r = 0.27–0.34, p < 0.05, n = 63 and r = 0.27–0.47, p < 0.05, n = 52, respectively) from the premonsoon period to the end of monsoon season (April–October).

More significant negative correlations were observed with almost all of the temperature indices, except for the monthly average minimum temperature, which showed no significant correlation with ∆¹³C chronology throughout the entire year. On a regional scale (CRU TS 3.24), monthly mean temperature exhibited negative correlations in March and May–October (r = –0.46 to –0.19, p < 0.01, n = 113); it exhibited the highest correlation in August, which was the month with the highest amount of rainfall. DTR and T_max exhibited similar correlations with ∆¹³C chronology from January to September, except in April (i.e. the hottest month of the year), in which T_max exhibited no significant correlation and DTR had the lowest correlation among 10 correlated months. The correlation coefficients between DTR and ∆¹³C chronology ranged from –0.19 to –0.47 (p < 0.01, n = 113), and those between T_max and ∆¹³C chronology ranged from –0.26 to –0.58 (p < 0.01, n = 113). On a local scale (MHS meteorological station), tree-ring ∆¹³C values were significantly negatively correlated (p < 0.01) with monthly mean temperatures from April to October (i.e. r ranged from –0.27 to –0.54). Tree-ring ∆¹³C values were significantly negatively correlated (p < 0.01) with T_max values from April to June and August, as their correlation coefficients ranged from –0.30 to –0.47; they were also significantly negatively correlated (p < 0.01) with DTR values from January to June and August, as their correlation coefficients ranged from –0.24 to –0.45. Although ∆¹³C chronology was still significantly negatively correlated with the MHS meteorological DTR and T_max data in August (the month with the highest precipitation in the year) in the same way that it was on a regional scale, ∆¹³C chronology exhibited a better negative correlation with DTR and T_max during the premonsoon season (i.e. April to June, with the highest correlation on May) than it exhibited with DTR and T_max on a regional scale, which performed better throughout the whole monsoon season.

After combining data from each climatic parameter over several months to determine the local climate response (MHS meteorological station), ∆¹³C chronology was negatively correlated with the average mean temperature from the premonsoon period to the end of the monsoon period, i.e. from May to October (r = –0.60, p < 0.001, n = 63) (Fig. 4). For the regional climate response (CRU 3.24), ∆¹³C chronology was negatively correlated with the average maximum temperature during the peak of the monsoon period, i.e. from June to September (r = –0.62, p < 0.001, n = 113) (Fig. 4). Previous studies of cambium activity in two species of Thai pines (Pinus kesiya and Pinus merkusii) in the northern region of Thailand (Pumijumnong and Wanyaphet, 2006) revealed that the growth of pines started in May and ended in the early dry season (December to January). However, depending on the soil moisture level in the area, the cambium in some trees may still be active through the late dry season (March to April). Consequently, the δ¹³C values in tree rings should reflect the temperature during their growth period.

Although the Pinus merkusii tree-ring δ¹³C data in this study captured the climatic signal of the regional T_max and local T_mean data from May to October, which are related to the growth period, it is still unclear if the cambium activity of the trees in this study ceases in October, which means that the whole ring may potentially reflect the temperature during the time of tree growth, or if the cambium remains dormant throughout the dry season, which means that the whole ring reflects the temperature during the time that the early wood was formed. In the latter case, the δ¹³C value that is preserved in the whole ring did not exhibit significant correlation with the climate signal from November to April, which may be explained by the different proportions of wood in the two different periods (i.e. the monsoon period and dry period). Because the number of cells in the cambium zone is higher from May to October (i.e. in the monsoon period of the current year) than it is from November to April (i.e. during the following dry period) (Pumijumnong and Wanyaphet, 2006), it seems likely that the δ¹³C values preserved in the higher proportion of wood dominated the δ¹³C values preserved in the lower proportion of wood.

In general, the maximum temperature, especially during the monsoon season, is more likely to increase with decreasing rainfall because in cloudless conditions or clear weather, daytime temperatures tend to be higher as they receive more direct radiation from the sun. Therefore, the decadal warm period in the T_max series might reflect drought conditions. To determine whether this hypothesis is consistent with the MHS climatic data, we performed the Pearson’s correlation between the actual June–September maximum temperature (CRU TS 3.24), actual June–September average precipitation (CRU TS 3.24) and June–September average relative humidity (MHS meteorological station). All of the climatological data were calculated as a smoothed 11-year moving average. The fact that T_max exhibited negative significance with both precipitation (r = –0.5, p < 0.01) and relative humidity (r = –0.75, p < 0.01) confirmed this hypothesis; in particular, its high correlation with relative humidity directly reflected drought conditions. Therefore, we decided to reconstruct June–September maximum temperature data based on the relationship between ∆¹³C chronology and the CRU TS 3.24 June–September T_max, which had the strongest relationship among other climate signals and an indirect relationship between maximum temperature and other climate signals.

3.3. Transfer function and reconstruction of June–September mean maximum temperature data

Based on the relationship between ∆¹³C chronology and the CRU TS 3.24 June–September T_max, the following simple linear regression model was developed:(3) $$T_{\text{max}\text{Jun}-\text{Sep}}=-0.739{\mathrm{\Delta }}^{13}\text{C}+42.156$$(3)

(N = 113, r = –0.62, R² = 0.384, R²_adj = 0.378, F = 69.178, p < 0.0001)

where r is the correlation coefficient between T_{max Jun–Sep} and ∆¹³C chronology. The explained variance of the reconstruction was 38.4% (37.8% after adjusting for the degree of freedom). A comparison of the reconstructed and actual regional June–September T_max is shown in Fig. 5. The split-cross calibration and verification were performed to test the stability and reliability of the regression equation (Table 2). The values of the reduction error (RE) and the coefficient of efficiency (CE) were positive, which indicated that T_{max Jun–Sep} model potentially produced a legitimate reconstruction. Based on the linear regression model (Equation (3(3) $$T_{\text{max}\text{Jun}-\text{Sep}}=-0.739{\mathrm{\Delta }}^{13}\text{C}+42.156$$(3) )), the 226-year mean maximum temperature from June–September was reconstructed to cover the period from 1788 to 2013 AD. The mean T_{max Jun–Sep} value was 30.09 °C, although it ranged from 29.52 to 30.94 °C, with a standard deviation of ±0.29 °C (Fig. 6a).

Fig. 5. Comparison of the reconstructed and actual June–September maximum temperature from 1901–2013.

Table 2. The split-cross calibration of the June–September maximum temperature reconstruction.

Calibration	r	Sign test	Verification	r	RE	CE	Sign test
1901–1957	0.70**	38+/19–*	1957–2013	0.64**	0.272	0.229	40+/16–**
1957–2013	0.64**	42+/15–**	1901–1957	0.70**	0.132	0.042	36+/20–*
1901–2013	0.62**	80+/33–**

* p < 0.05.

** p < 0.01.

Fig. 6. The June–September maximum temperature reconstruction (grey line) with a 0.1-Hz low-pass filter (red line). (a) The standard deviation (σ) was ±0.29. (b) Four shifting segments calculated using the Bernaola-Galvan Segmentation Algorithm.

4. Discussion

4.1. δ¹³C variability and its climatic implication

In the δ¹⁸O-based precipitation reconstruction in northwestern Thailand (Xu et al., 2015), the δ¹⁸O time series from the 4 Pinus merkusii trees exhibited less inter-tree variability than the ∆¹³C series in our study, which was expressed in both its greater values of Rbar (0.66–0.83) and EPS (0.88–0.94) obtained from the same amount of samples. This greater inter-tree variability in δ¹³C values may be related to the different mechanisms that control the δ¹⁸O and δ¹³C values in tree-rings (McCarroll and Loader, 2004; Leavitt, 2010). Tree-ring cellulose δ¹⁸O values obviously reflect the δ¹⁸O values of ambient precipitation, which are generally dominated by precipitation and temperature; meanwhile, tree-ring cellulose δ¹³C values are controlled by the balance between the stomata conductance rate and the photosynthesis rate. Stomatal conductance is influenced by the factors that control ambient water (i.e. precipitation, relative humidity, and soil moisture status, which are generally negative correlated with δ¹³C), and the photosynthesis rate, which is influenced by the factors that control the ambient temperature and light level (i.e. mean, minimum and maximum temperature and solar irradiation, which are generally positive correlated with δ¹³C).

In our results, the ∆¹³C values of Pinus merkusii showed weak correlations with hydroclimatic factors (which modify the stomata conductance rate) but exhibited strong correlations with temperature (which modifies the photosynthetic rate). Many previous studies have reported weak correlations between precipitation and δ¹³C, thus suggesting that tree-ring δ¹³C values are less sensitive to precipitation in non-arid areas where ambient water is not a main limiting factor (Warren et al., 2001; Porter et al., 2009; Bégin et al., 2015). The annual precipitation of 1228 mm in our study area supports this hypothesis; although 90% of the annual rainfall occurs from May to October, the average monthly relative humidity during the dry season (November to April) is 68%, which keeps this area far from experiencing moisture stress during the dry season.

4.2. Characteristics of the June–September maximum temperature reconstruction

To investigate the physiology of the non-stationary T_maxJun–Sep reconstruction under shifting climatic conditions, the Bernaola-Galvan Segmentation Algorithm (Bernaola-Galván et al., 2001) was run on the T_{max Jun–Sep} reconstruction (Fig. 6b). The following four lengths of segments with different local mean T_{max Jun–Sep} values were revealed to be significant at the 0.01 level: 1788–1829 (mean = 29.83 °C), 1830–1869 (mean = 30 °C), 1870–2001 (mean = 30.22 °C) and 2002–2013 (mean = 29.91 °C). During the first two segments, from the late 1700s to the mid-1800s (i.e. 1788 to 1869 AD, which are assumed to be the pre-Industrial Revolution period), the maximum temperature increased slightly at a rate of 0.03 °C/decade, but the mean maximum temperature was obviously lower than that of the third segment (1870–2001), which represents the post-Industrial Revolution period (29.91 °C vs. 30.17 °C). During the third segment (1870–2001), the maximum temperature pattern seemed to be constant compared to the changing rate (+0.004 °C/decade). However, the temperature fluctuations were rather high compared with that of the previous time period (σ ± 0.28 vs. σ ± 0.20). The short fourth segment, which occurred from 2002 to 2013, showed a deceasing trend at a rate of 0.12 °C/decade.

We defined a warmer year as >mean +σ (30.38 °C) and a cooler year as <mean – σ (29.8). We defined 34 warmer years and 40 cooler years, which accounted for 14.91 and 17.70% of the entire reconstruction period, respectively. Interestingly, all of the ten warmest years appeared after 1900 (from highest to lowest temperatures, these years include 1950, 1949, 1948, 1947, 1945, 1946, 1987, 1998, 1989 and 1979), in contrast to the ten coolest years, which were scattered throughout the entire series (from lowest to highest temperatures, these years include 1938, 1807, 1940, 2006, 1821, 1817, 1789, 1936, 1797 and 1939). Consequently, the Mae Hong Son tree-ring δ¹³C–base maximum temperature data potentially captured the effects of twentieth century global warming, which may not exhibit obvious increases in maximum temperature trends but may instead reflect the increased severity of temperature fluctuations.

To observe decadal maximum temperature patterns, 11-year moving averages were obtained in order to emphasize low-frequency climate signals. Four cool and three warm periods were defined by average temperatures that were lower and higher than the 228-year mean temperature, respectively. The cool periods were 1788–1871, 1931–1939, 1956–1974 and 2004–2013. The three warm periods were 1872–1930, 1940–1955 and 1975–2003. It should be noted that all three of the warm periods appeared after 1850 AD.

The mean temperature from the first cool decades (1788–1829) in the T_max reconstruction is the lowest among all four of the cool periods, with a mean maximum temperature of 29.82 °C. This condition may result from a negative climate forcing phase. Negative solar forcing and volcanic forcing (which is also known as volcanic-solar downturn) during 1791–1820 has been reported in volcanic forcing reconstructions based on ice core index analyses (Gao et al., 2008; Crowley and Unterman, 2013) and solar irradiance reconstruction (Shapiro et al., 2011; Vieira et al., 2011). Two negative climate forcing events coincidentally occurred during the same period, namely, the Dalton minimum (which featured low solar activity due to a low sunspot count) from 1800–1820 AD (Shapiro et al., 2011) and the eruptions of two volcanoes in 1809 and 1815 (Gao et al., 2008); these events caused a significant decrease in global temperature.

The three warm periods in our monsoon maximum temperature reconstruction, especially those that occurred during the late nineteenth century, coincided with drought periods that have been observed in many tree-ring-based studies. Pumijumnong and Eckstein (2010) reconstructed premonsoon temperatures based on Pinus merkusii and Pinus kesiya tree-ring chronology along northwestern Thailand that reflected drought periods during 1880–1910, 1950–1965 and 1980–1990. Buckley, Palakit, et al. (2007) used teak chronology from Mae Hong Son as a proxy for May–December drought conditions. Mid-eighteenth- and late-nineteenth-century droughts were reflected in this chronology. Additionally, the March–May PDSI reconstruction over northern Vietnam based on Forkienia hodginsii tree-ring widths revealed two prominent droughts that occurred in the mid-eighteenth and late nineteenth century (Sano et al., 2008).

4.3. Teleconnection with climate driving force

The coupled ocean–atmosphere system controls the characteristics of the monsoon climate over Thailand. This region is compressed on two opposite sides by the Pacific Ocean and the Indian Ocean, which produce mass driving mechanisms that have broad-scale influence. Generally, monsoon rainfall in Thailand is a result of the interplay between the south-west monsoon, the sub-region of the Indian summer monsoon and the movement of the Inter-Tropical Convergence Zone (ITZC), which brings a warm-moisture air mass from the Indian Ocean towards Thailand from May to October. The Indian Ocean Dipole (IDO), which is an irregular warm–cool oscillation sea surface temperature in the western Indian Ocean, represents a Dipole mode index (DMI) (Saji et al., 1999), which is a sea surface temperature index that can be used to calculate the difference between the Western Tropical Indian Ocean index (WTIO) and the southeastern Tropical Indian Ocean index (SEITIO); this was used to investigate the influence of the Indian monsoon on the temperature over northwestern Thailand. Although the precipitation over northwestern Thailand is mainly affected by the Indian summer monsoon, the values of the DMI, WTIO and SEITO indexes showed no significant correlation with the Tmax reconstruction.

The effects of ENSO on the inter-annual and inter-decadal variability of temperature over Thailand have been reported, as El Niño tends to exhibit higher temperatures than normal (with the opposite occurring during La Niña events) (Limsakul and Goes, 2008). Similar to the results of our T_max reconstruction, high peak temperature values were observed in 1919, 1941, 1958, 1972, 1982, 1987, 1995, 1998 and 2009 and they coincide with previous winter and/or current El Niño years. Conversely, no obvious T_max reconstruction showed low temperature peaks during La Niña years. Thus, the centre Pacific El Niño and the eastern Pacific El Niño exhibit discrete influences on monsoon precipitation (Weng et al., 2009). Therefore, the HadlSST Nino 3 index (1870–2013) was selected to represent the eastern Pacific El Niño, and the HadlSST Nino 4 (1870–2013) was selected to represent the centre Pacific El Niño (Weng et al., 2009; Yeh et al., 2009; Xu et al., 2015). The Pearson’s correlation between our maximum temperature reconstruction and the two Nino indexes exhibited weak correlation, but it was significant during the beginning of the monsoon period (May–June), in which r_{Nino 3 vs. Tmax} was 0.21 and p < 0.05 and r_{Nino 4 vs. Tmax} was 0.19 and p < 0.05. To investigate changes in the long-term relationship between T_max reconstruction and the Nino index, a 21-year running correlation was performed. In Fig. 7, the 21-year running correlation between T_max reconstruction and the 2 Nino indexes exhibited the same positive unstable correlation trends. Low correlations appeared during 1899–1958 and 1984–2013, while higher correlations appeared during 1788–1898 and 1959–1983. The fact that the relationship between Nino 3 SST and Indian summer monsoon rainfall has weakened since 1980 (Kumar et al., 1999) may be a reason for the decreasing correlation trend between T_max reconstruction and the two Nino indexes observed after 1984 in our study. The south-eastward shift of the descending limb of Walker Circulation, which broke down the relationship between ENSO and Indian monsoon, had significantly different responses to precipitation between the northwestern and central regions of Thailand. Decreasing correlations after 1980 were also observed in the relationship between the tree ring δ¹⁸O-based July–October precipitation reconstruction in northwestern Thailand and ENSO (Xu et al., 2015), but ENSO has since shifted its influence to overcome the precipitation over central Thailand (Singhrattna et al., 2005).

Fig. 7. (a) 21-year running correlation between June–September T_max reconstruction and May–June Nino 3 and Nino 4; (b) May–June sea surface temperature of the PDO index.

The low–high correlation between T_max reconstruction and the two Nino indexes appears to reflect an oscillation pattern, which raises more questions about whether there are any possible phenomena in the Pacific Ocean that could have produced this pattern. This was investigated using the warm–cool phase (positive–negative phase) of the Pacific Decadal Oscillation (PDO) index. Although the May–June PDOersst index (https://www.ncdc.noaa.gov/teleconnections/pdo/) from 1880–2013 did not show significant correlations with the T_max reconstruction, a 21-year running correlation between T_max reconstruction and the PDO index tended to exhibit a positive correlation during periods in which the T_max reconstruction and the Nino index exhibited high correlations, and they exhibited a negative correlation during periods in which the T_max reconstruction and Nino index exhibited low correlations (Fig. 7a). Interestingly, the positive phase of PDO seems to decrease the relationship between T_max over northern Thailand and ENSO (and vice versa) (Fig. 7b). The effects of PDO and ENSO on each other have been explained in terms of their magnitudes. Gershunov and Barnett (1998) concluded that highly positive PDO values may lead to stronger El Niño events and, conversely, that stronger La Niña events tended to occur during periods of highly negative PDO values. However, it still is too early to conclude that PDO exerts an influence on T_max in northwestern Thailand or even physically controls ENSO; thus, further investigation is needed.

5. Conclusions

The stable carbon isotopes of the Merkus pines in northwestern Thailand significantly captured temperature signals over the core-monsoonal season. Thus, a 226-year June–September maximum temperature record extending back to 1788 AD was reconstructed based on ∆¹³C chronology. The analysis of climatic shifting conditions revealed a slightly increasing trend in maximum temperature before 1869 AD, followed by increasing temperature fluctuations. In addition, warm period from 1872–1930 observed in this reconstruction were consistent with the drought periods that have been observed during the late nineteenth century in tree-ring-based studies in Thailand and Vietnam. Unstable relationships were detected between the May–June Nino 3 and Nino 4 indexes and our reconstruction; in particular, a weaker relationship after 1980 may be related to a shift in the descending limb of the Walker Circulation. Furthermore, we discovered that this unstable relationship coincided with PDO; specifically, negative phases of PDO tend to increase this relationship and vice versa.

Disclosure statement

No potential conflict of interest was reported by the authors.

Funding

This study was jointly supported by Grants from the CAS Key Research Program of Frontier Sciences [grant number QYZDJ–SSW–DQC021], [grant number NSFC41630531], [grant number XDPB05], [grant number GJHZ1777]; the Key Project of IEECAS and the SKLLQG.