A Possible Volcanic Signal in NCEP/NCAR Reanalysis Data

In this study we used the multi-channel singular spectrum analysis (MSSA) to extract climate signals from National Centers for Environmen­ tal Prediction I National Center for Atmospheric Research (NCEP/NCAR) reanalysis monthly mean temperature data between 1950 and 1998. The signal associated with three major volcanic eruptions during this period was clearly shown in the sixth principal component (PC) of the MSSA. The time evolution of this PC revealed both interannual and interdecadal scale variations suggesting that major volcanic eruptions not only have an effect on short-term climate variability but may also be associated with interdecadal variability in the atmospheric temperature field. The associ­ ated anomaly patterns of this PC showed general cooling in the lower at­ mosphere and warming in the upper atmosphere. Furthermore its time lag maps show anomaly patterns in the equatorial Pacific similar to those of El Nino-Southern Oscillation (ENSO) warm events. These results suggest that a major volcanic eruption might enhance the strength of an ENSO warm event. Nevertheless, since these eruptions occurred rarely and since each one of them roughly coincided with an ENSO warm event, we are unable to assert that there is a relationship from observed data alone. Despite the irregular occurrence of major volcanic eruptions in history, this study demonstrated that MSSA has the ability to extract temporally irregular but meaningful signal from observed data. Hence MSSA is a very valuable tool for climate study. (


INTRODUCTION
Major explosive volcanic eruptions inject massive amounts of dust and gases into the 1 1nstitute of Atmospheric Physics, National Central University, Chung-Li, Taiwan atmosphere.These volcanic aerosols can reduce incoming solar radiation and cool the earth's surface and troposphere.Furthermore they can remain in the stratosphere for 2-3 years (Mass and Portman 1988).Therefore major volcanic activity is believed to have a significant impact on global as well as regional climate variability.However, the rare occurrence of major explo sive eruptions and the presence of other interannual variability in climate records, especially ENSO variability, render detection of the volcanic signal a complex issue.Numerous studies have attempted to extract such a signal from observational data.Angell (1988) was the first to address the problem that the volcanic signal may be masked by ENSO variability, and used lag correlation between tropic sea surface temperature and tropospheric temperature to subtract the ENSO signal from the tropospheric temperature time series.His results indicated that major volcanic eruptions cooled the northern hemispheric continental surface by about 0.3°C.Subsequently, Mass and Portman (1989) used the superposed epoch method of analysis (Composite method) to examine annual and monthly variations in four factors after major volcanic eruptions during the past century; these were surface air temperature, sea surface temperature, surface pressure, and precipitation.They used various data (such as the Volcanic Explosive Index and the Dust Veil Index) to determine if an eruption was a major one, but the main selection criterion was whether the eruption created a significant stratospheric dust veil.They established that there were nine major volcanic events between 1883 and 1983.They found that after the very largest of these eruptions (in terms of stratospheric dust cloud), there was a modest cooling in zonal mean surface air temperature.They also used regression coef ficients derived from tropical SST to remove the ENSO signal from data and found an en hanced volcanic signal in air temperature composited for the five largest eruptions.Roback andMao (1992, 1995) extended Mass and Portman's analysis to regional and seasonal pat terns with a global surface temperature dataset.They used regression coefficients derived from the Southern Oscillation index (SOI) to remove the ENSO signal from .surface tempera ture anomalies at each grid point.They found that warming occurs over Eurasia and North America, while cooling occurs over northern Africa, the Middle East, southern Asia, and Australia.More recently, Portman and Gutzler (1996) also extended Mass and Portman's analysis to detect volcanic signals in U.S. surface climate records.Instead of using regression relations to remove the ENSO signal, they separately composited major explosive volcanic eruptions, ENSO warm events, and ENSO cold events occurring since 1900.Then they sub tracted the composited ENSO-related anomalies from the anomalies of temperature and pre cipitation associated with the volcanic eruptions.After removing the ENSO-related signal they were able to detect volcanic signals in monthly records; these showed that warming oc curs east of the Continental Divide and cooling occurs west of Continental Divide and the southern United States.
The rare occurrence of major explosive eruptions meant that past detection of volcanic signals relied heavily on the superposed epoch method of analysis.Moreover, the frequent coexistence of ENSO events with major eruptions also forced these studies to assume that these two signals are independent to allow the removal of the ENSO signal by simply subtract ing it from data.This kind of simple subtraction procedure may artificially suppress a volca nic signal in strong ENSO response regions.Cons�quently, results from the above studies may not present the complete picture of the volcanic signal in the atmosphere.Hence a more objective and less restrictive method is needed to extract a volcanic signal from observed data to reveal its true identity.
According to Mass and Portman (1989), a volcanic signal should include the following characteristics: first, it should first appear after an eruption; second, it should be observed within 2-3 years after an eruption; third, it should maintain the same sign and approximate magnitude for several months; and fourth, a multieruption composite should be similar to the individual component eruptions in magnitude and sign.These characteristics suggest that a volcanic signal have an extensive spatial-temporal scale.To objectively extract the volcanic signal from observed data, the widely used multi-channel singular spectrum analysis (MSSA; Vautard et al. 1992;Plaut and Vautard 1994;Robertson et al. 1995;and Fraedrich et al. 1997) has appeared in recent years to be a valuable tool.This is because it takes lag covariance among time series into consideration, and therefore can effectively extract oscillation patterns and can differentiate variability associated with different time scales in data.However, for MSSA to reliably extract its signal, the rare occurrence of major explosive eruptions also requires data to have extensive spatial-temporal coverage.Fortunately, the recently available National Centers for Environmental Prediction I National Center for Atmospheric Research's (NCEP/NCAR) reanalysis dataset (Kalnay et al. 1996) has the necessary consistent quality as well as the extensive spatial-temporal coverage of many atmospheric variables that is necessary.
This paper is arranged as follows.In section two, we briefly describe the data and analy sis methods used in this study.In section three, we present results of the MSSA analysis.Finally, in section four, we summarize the major findings and discuss some of the implications.

DATA AND ANALYSIS MET HOD
The original data is the NCEP/NCAR's reanalysis monthly mean temperature field be tween 1950 and 1998.In the analysis, the horizontal domain and resolution are 60°S to 60°N and 5° longitude by 5° latitude, respectively.In the vertical, twelve levels of data are included in the analysis; these are 1000, 925,850,700,600,500,400,300,250,200,150 and 100 hPa respectively.
The principal analysis tool in this study is the MSSA method.In short, the method is an extension of the principal component analysis (PCA; Preisendorfer and Mobley 1988) of a time series of spatial vectors, including, in addition, an analysis of temporal structure.This method first expands a vector of K components that evolves with time into a new vector of K x W components by sliding a time window of length W onto the original vector, then performing a PCA on this new where m is the length of time series and K is the total channel number, then by sliding a time window of length W onto the original vector, the expanded data matrix takes the form After Z has been constructed, we can perform a conventional PCA analysis, or alternatively, we can directly apply singular value decomposition (SVD; Golub and Van Loan 1989) to Z to obtain associate eigenvalues, principal components (PCs) and empirical orthogonal functions (EOFs).Vautard et al. (1992) noted that because the MSSA takes lag covariance among time series into consideration, it can effectively differentiate the variability associated with differ ent time scales in data.Furthermore, they showed that a pair of nearly identical eigenvalues indicates the existence of an oscillation pattern in data.Therefore the MSSA has the ability to extract quasi-periodic behavior from the data.
Before the MSSA, we first removed the climate mean annual cycle from data at each grid point and normalized the resultant anomalies using their standard deviations.Then, at each level, the normalized data was subjected to a PCA; only the first 60 PCs was retained to compress the data.Next, we gathered all of these PCs to form a multi-variable matrix and then performed another PCA on this matrix.Note that this step is essentially equivalent to the multi-variable empirical orthogonal function analysis (MV-EOF).Wang (1992) suggested that, by considering covariance among variables, the MV-EOF was therefore able to enhance the signal-to-noise ratio and efficiently extract recurrent patterns in the data.Because it is very unlikely that dominant atmospheric variability can occur only in a single level, we can expect that MV-EOF analysis of the twelve levels in this study will effectively enhance the true signals of the data.To further compress the data, we only used the resultant first 30 PCs with a 30-month window length to perform the MSSA analysis.
To test the statistical significance of MSSA' s results, we adopted the Monte-Carlo MSSA from Zhang (1998), which originated from Allen and Smith (1996), to determine if the EOFs of MSSA are significantly different from what could be obtained from pure red noise.This is because, as Allen and Smith (1996) noted, temporal autocorrelation may cause high-variance low-frequency modes of MSSA that entirely due to noise to be identified as a signal.There fore a MSSA result should be tested against a red noise null hypothesis to evaluate its statisti cal significance.In this study we used a red noise model, which consists of 30 independent AR(l) processes, each with the same variance and first-order autocorrelation as those of the first 30 PCs of the MV-EOF, to generate 200 segments of noise data with the same sample size as the observation.The noise data was then projected onto the EOFs of MSSA to obtain the amount of variance explained by each of these EOFs.If an eigenvalue of the MSSA is greater than the second (the fifth) largest variance of the 200 noise data associated with the same EOF, then we would have 99% (95%) confidence that this mode is significantly different from red noISe.
One of the advantages of MSSA is that its EOFs contain both space and time information (which are thus commonly called ST-EOFs).Therefore we can not only examine each PC's Zs be the sth ST-EOF and Ykqs the associated lag q spatial pattern at the k1h level, where Then, according to (2), y kq s can be derived as In the following section, all the spatial patterns were constructed following the procedure from (3) to (7).

3.RESULTS
Figure 1 is the percentage of the total normalized NCEP temperature data explained by the first 10 modes of the MSSA.Also included are the 99% and 95% confidence intervals from the Monte Carlo MSSA of pure red noise.Note that the fourth, sixth, seventh, ninth and tenth modes have passed the 1 % significance level.Therefore we have very high confidence that these modes are signals of the NCEP data.If a 5 % significance level is used, then the second and the third modes can also be considered as signals.This is because, in MSSA, two nearly identical eigenvalues implies that these two modes form a pair of oscillations.Therefore, even though mode 3 does not pass the 5% significance level, because its eigenvalue (6.6%) is very close to that of mode two (6.9%), it should also be considered a signal.Note from this figure that although the first mode explained more than 23% of the total variance, it was not significantly different from red noise at 1%or5% level.Thus a high-variance mode in MSSA does not necessarily mean that there is high confidence that the mode is a signal.
Figure 2 (a) to (f) shows respectively the first to the sixth MSSA PCs of the NCEP data.
Although mode one did not differ significantly from red noise, its high variance still deserves attention.Figure 2(a) shows how the first PC for the period before 1975 is opposite in sign to that for the period after 1980, with a sharp transition in the late 1970s.This suggests that the  first mode may be a multi-decadal mode or may correspond to the wanning trend that began in the late 1970s (Graham 1995).The second and the third PCs correspond very well to the occurrence of major ENSO events between 1950 and 1996 (plus signs in the figure mark the initial months of ENSO warm events, as defined in Trenberth, 1997).Thus they constitute a pair of oscillations to describe the 3-5 year evolution of ENSO variability.As for the fourth and fifth PCs, they appear to be a mixture of signals from a 2-3 years ENSO variability and a low frequency mode.
The. sixth PC in Fig. 2(f) has the most interesting feature.It has three sharp peaks that correspond very well to the dates of three major volcanic eruptions.These were (Mass and Portman 1989) Agung (March 1963), El Chichon (April 1982) and Pinatubo (June 1991) and are each marked with a cross signs in the figure.Furthermore, each of these sharp peaks is followed by a long oscillatory downward period.This behavior matches quite well with gen eral characteristics of a volcanic signal, as described in the introduction.To test the robustness of the above features, we have repeated the analysis with both different window lengths and different analysis period lengths.For window lengths from 6 month up to 60 months, results are very similar.Because the volcanic effect is strongest during the first year after the eruption, window length larger than 60 months will apply too much filtering on high frequency variability.This prevents the volcanic signal from being correctly shown in the leading modes of the analysis.Considering different analysis periods; there is 19 years between the Agung and El Chich6n eruptions, while only 9 years between El Chich6n and Pinatubo eruptions.It can be seen that the shorter the period that includes only the latest two eruptions, the stronger the volcanic signal.For instance, Fang ( 1998) also analyzed the same NCEP/NCAR data for the period 1982 to 1994 and found a possible volcanic signal presented in the fourth mode of his Jan50 Jan60 Jan70 Jan80 Jan90 Jan50 Jan60 Jan70 JanBO Jan90 (e) ...----------� -0.1 Jan50 Jan60 Jan70 Jan80 Jan90 Jan50 Jan60 Jan70 Jan80 Jan90 Fig. 2. (a) to (t) show the first to the sixth PC respectively.Plus signs denote the beginning dates of ENSO warm events as defined in Trenberth (1997).
Crosses are the eruption months of three major volcanic eruptions in the period.Since the 30-month window length was used, these PCs were shorter than original data.Therefore they ended at July 1996.

901
MSSA.Nonetheless, since El Chich6n and Pinatubo eruptions happened to coincide with ENSO warm events, without additional evidence from the Agung eruption, we cannot com fortably identify the fourth mode in Fang (1998) as a true volcanic signal.On the other hand, this study used longer period.But the occurrence of only three events in 50 years also inevita bly leads to a weakening of the signal.This dilemma is the main difficulty in extracting clear volcanic signal from observed data.Nevertheless, the almost exact match of the major peaks with the dates of three major eruptions in the analysis period suggests that the sixth PC is indeed the volcanic signal.Furthermore, if signals can be detected in leading modes of MSSA for a mere three events in 50 years, this suggests that major volcanic eruptions may have a significant impact on climate.
To examine the regional impact of a volcanic signal, we constructed the spatial patterns of the sixth PC from its associate EOF. Figure 3 shows its zero-month lag spatial patterns at 1000, 500 and 150 hPa.In the lower troposphere, we can see that most regions are occupied by negative anomalies.There are major negative anomaly centers near the Philippines, on the  and Mao (1995), where the ENSO signal was subtracted from the data, our results show major differences in the tropical regions, where there was a strong ENSO signal.Be cause the three volcanoes studied are all located in tropical latitudes around the Pacific (Agung, 8°S; El Chichon, 17°N; Pinatubo, 15°N), the existence of a strong volcanic signal in the equa torial Pacific is not unexpected.In fact, a general circulation model simulation (Robock and Liu 1994) also showed a strong volcanic signal in the tropics.Thus these results suggest that the removal of the ENSO signal by simple.subtractionmay artificially suppress a volcanic signal in strong ENSO response regions.

Roback
In the mid-troposphere, most of the cooling occurs in tropical latitudes with centers near the coast of South America at around 15°S, Central America, the north Pacific near dateline and equatorial Africa.On the other hand, warming mainly occurs in mid-latitudes with cen ters near southern Australia, the tip of South America, the Aleutian Islands, northern China and the British Isles.Compared with patterns at 1000 hPa, it is interesting to note that the anomaly centers near the western coast of South America have changed in sign from positive to negative, while negative anomaly centers near the Philippines and the Middle East have almost disappeared.As for the upper troposphere and lower stratosphere, most areas experi enced warming or no change.The cooing area was mainly confined to the equatorial western Pacific, Australia and the American continent.This is consistent with the effect of volcanic dust in the upper atmosphere.
Figure 4 is the same as Fig. 3 except for the six-month lag spatial patterns.In most regions, the anomaly patterns are similar to those in Fig. 3.However, in the equatorial Pacific at 1000 hPa, there are some obvious differences.Here, the warm anomaly now occupies the whole central and eastern Pacific, while the cold anomaly extends from the maritime continent northeastward and southeastward toward the subtropical central Pacific.These features are similar to those of the ENSO mode (Rasmusson and Carpenter 1982;Harrison and Larkin 1998).For the sake of comparison, Figure 5 shows the zero-month lag spatial patterns of the second PC (with sign reversal) that corresponds to the spatial distribution of ENSO warm events.Note that because PC 2 and PC 3 form a pair of oscillation modes to describe the main evolution of the ENSO, the phases of their ST-EOFs are in quadrature, but their spatial pat terns are quite similar.Therefore either one of them is adequate to show the spatial distribu tion of the ENSO.We can see, in the equatorial Pacific at 1000 hPa, that they indeed share some common features.Nevertheless, the warm anomaly centers are located in different regions; i.e., the maximum warming regions are located near the coast of South America and near the Hawaiian Islands in Fig. 4, but near 140°W at the Equator in Fig. 5. Furthermore, in regions other than the equatorial Pacific at 1000 hPa as well as in 500 and 150 hPa, they are quite different.Thus although the six-month lag spatial patterns of mode 6 in the equatorial Pacific at 1000 hPa bear some similarities with the ENSO mode, the differences in other regions suggest that these two modes are distinct.
Figure 6 shows the 24-month lag spatial patterns of the sixth PC.We can see that the anomaly patterns in the equatorial Pacific are still clearly identifiable on the map.The ENSO mode, on the other hand, has already changed to a different pattern on the 12-month lag map (not shown).These features suggest that the spatial patterns of the sixth PC persist longer than a typical ENSO oscillation.Therefore it appears that this mode is indeed not associated with  Fig. 6.As in Fig. 3, except for the 24-month lag spatial patterns for the sixth PC.
the tropical Pacific suggest that the sixth PC of the MSSA may indeed be the volcanic signal.
Its time evolution suggests that major volcanic eruptions may not only cause dramatic tem perature change on an annual time scale but may also induce change in the decadal time scale.
Nevertheless, the use of a 30-month window length prevents us from ruling out the possibility that the sixth mode may be contaminated by interdecadal ENSO-like variability� Therefore the influence of major volcanic eruptions on a decadal scale is still an open question.After major t ropical volcanic eruptions, time-lag spatial maps of mode 6 show that the surface temperature field in the tropical Pacific is similar to that of ENSO warm events.There fore it is likely after these major tropical volcanic eruptions, enhanced cooling in the tropical western Pacific may cause significant changes in Hadley and Walker .circulations to induce a persistent surface temperature response similar to those of ENSO warm events.These results sµggest that major volcanic eruptions may enhance the strength of an ENSO warm event.Furthermore it appears that the removal of the ENSO signal by simple subtraction may indeed artificially suppress the volcanic signal in strong ENSO response regions.Nevertheless, be cause each of these three eruptions more or less coincided with an ENSO warm event, it is almost impossible to ascertain this effect from observed data alone.Thus a detailed general circulation model simulation is needed to verify this possibility.
From this study, a rule of thumb can be made regarding the choice of a proper window length to allow MSSA to successfully extract an irregular signal from observed data.Note that the spectrum resolution of an MSSA strongly depends on the window length, i.e., a larger window length corresponds to a high spectrum resolution for the MSSA decomposition and vice versa.If the window length chosen too large compared to the dominant variation time scale of an irregular signal, then the irregular nature of the signal will likely be ranked behind more regular signals and hence go unnoticed.On the other hand, if the window length chosen too small, then the low spectrum resolution is likely to cause several signals of various time scales to mix together to allow clear identification of the signal.Thus, i. t appears that a proper window length should be based on the dominant variation time scale of the signal that one is interested in.This study demonstrated that even though major volcanic eruptions have oc curred irregularly in history, MSSA has the ability to extract temporally irregular but mean ingful signals from observed data.Hence MSSA is a very valuable tool for climate study.

j=l
associated spatial patterns, but also their time evolution from the associated ST-EOF.To extract the spatial-temporal structures of ST-EOFs, we need to construct spatial patterns of MV-EOF first.Let ejk be the jth EOF at the k'h level and br the r1h EOF of the MV-EOF, where N=l800 is the total number of grids and p = 60 X 12 is the total number of PCs used in MV EOF.Then the spatial pattern at the k'h level of the rth EOF of the MV-EOF, xkr, Once xkr has been derived, we can proceed to construct spatial patterns of the ST-EOF.Let

Fig. I .
Fig. I.The percentage of the total normalized NCEP temperature data explained by the first 10 modes of MSSA.Also included are the 99% and 95% confidence intervals from the Monte Carlo MSSA of pure red noise.
The percentage of total variance explained by the first 1 O modes of M SSA