Measurement of the stochasticity of low-latitude geomagnetic temporal variations

Ground magnetometer measurements of total magnetic field strength from 6 stations at low latitudes were analyzed using power spectrum and Hurst range scaling techniques. The Hurst exponents for most of these time-series were near 0.5, which indicates stochasticity, with the highest latitude stations exhibiting some persistence with Hurst exponents greater than 0.6. Although no definite correlations are evident, the relative increase of the Hurst exponent with latitude suggests the possibility that the underlying dynamics of the magnetosphere change with latitude. This result may help quantify the dynamics of the inner magnetosphere itself without the direct presence of the solar wind driver.


Introduction
Strong nonlinear coupling between the solar wind and the earth's magnetosphere results in many dramatic disturbances in the near-earth space environment, such as the dynamic magnetic and auroral signatures as well as magnetotail plasma signatures associated with magnetic storms and magnetospheric substorms. The strongest coupling between the solar wind and the magnetosphere occurs near the magnetopause, close to magnetic field lines that map to the high latitude ionosphere. The irregular nature of high-latitude disturbances, which typically occur above 60 • geomagnetic latitude, are clearly manifested in the auroral electrojet indices AE and AL. Studies of these indices have suggested that the magnetosphere behaves as a self-organized system with a small number of degrees of freedom (Vassiliadis et al., 1990;Sharma et al., 1993) although there are questions as to whether the magnetosphere itself is a self-organized system or whether it simply reflects the self-organized state of the solar wind and interplanetary magnetic field (IMF) to which it is strongly coupled (Price and Newman, 2001).
Correspondence to: J. A. Wanliss (james.wanliss@erau.edu) If the magnetosphere is self-organized, it would imply that only a small number of coupled, nonlinear ordinary differential equations are required to describe its dynamic behaviour. Indeed, the number of degrees of freedom that describe the high-latitude phenomena was found to have an average value of about 3.3 (Vassiliadis et al., 1990;Roberts, 1991;Shan et al., 1991;Sharma et al., 1993), although these results have been challenged by other studies Price, 1992, 1993). Such nonlinear approaches have led to the development and improvement of various dynamic models of substorms (Ohtani et al., 1995;Baker et al., 1997Baker et al., , 2000Takalo et al., 1999;Horton et al., 2001) Whereas the studies cited above have focused primarily on magnetic variations due to high-latitude current systems, in this paper we consider low-latitude magnetic variations ( ∼ 35 • − 40 • ). The high-latitude magnetic field lines are connected to magnetospheric regions that map very closely to the solar wind but the low-latitude magnetic field lines are connected directly to the inner magnetosphere. Since the field lines at low latitudes (L ∼ 1.5 − 1.8) are almost dipolar, they are not as strongly influenced by the interplanetary medium as the high-latitude regions where chaotic signatures might simply reflect similar solar wind conditions. Thus, examination of low-latitude ground magnetometer signals can provide clues as to whether the magnetosphere is inherently self-organized.
At lower latitudes, the dominant magnetic variations are due to the two diurnal solar quiet (Sq) large-scale current systems with foci at about 30 • magnetic latitude in the ionosphere and peak current densities near local noon. Superimposed on these diurnal signals are higher frequency variations from magnetospheric sources. Coupling of the lowlatitude regions with the magnetosphere is achieved along magnetic field lines that map to the inner magnetosphere, and through the variation of the ring current at distances of about 3-5 earth radii (R E ), especially during magnetic storms. The references quoted above have demonstrated a clear nonlinear behaviour of high-latitude time-series, but the latitudinal extent and variability of such behaviour is unknown. We are, therefore, interested in determining whether the quantitative nonlinear dynamics of the high-latitude regions extends to lower latitudes; for example, to determine the fractal dimension pertinent to low-latitudes and whether this is similar to that found in the high-latitude studies. This paper represents a preliminary effort to address the issue. A global index that characterizes the low-latitude magnetic variation, is D st . However, since D st is computed only at hourly intervals, it was deemed advantageous to use actual magnetometer records that provide a higher time resolution in order to investigate the fractal properties present at low-latitudes. Furthermore, the individual magnetometer time-series give localized estimates of the fractal parameters rather than the global output that is sampled by geomagnetic indices. We have analyzed individual magnetometer records for five days in January 1993 during which magnetospheric activity indicates no magnetic storms although several small substorms are observed at high-latitudes. We find evidence that the dynamics of low-latitude regions of the magnetosphere, sampled by the magnetometer stations in the study, are primarily stochastic, although two stations exhibit signals that are not inconsistent with self-organized criticality but with a lower level of complexity than for the dynamics governing high-latitudes; that is, the calculated fractal dimension is lower than that found in the high-latitude studies cited above. Although the degree of coupling with the solar wind and the IMF is not clear, our results do suggest that the low-latitude inner magnetosphere behaves in a fundamentally different way to the high-latitude regions.

Data
We utilize data from an Anglo American Corporation experiment, which recorded geomagnetic temporal variations from an array of six magnetometer stations spread in latitude over South Africa (Wanliss, 1995), shown in Fig. 1. The stations were at Ellisras (Ellis), Bronkhorstspruit (Bronk), Lanseria (Lans), Vryburg (Vry), Boshof (Bos) and Hermanus (Herm), whose geographic locations, corrected geomagnetic latitudes and L-shell positions are listed in Table 1. The motivation for the experiment was to provide a rigorous understanding of the background magnetic field in the region, to be used in the interpretation of aeromagnetic surveys. A useful byprod- uct was the high time resolution data used in the present study. Simultaneous data from the stations were acquired during 13-18 January 1993. The magnetometer instrument sensitivity is 0.1 nT with a sampling interval of 10 s. The total ambient field in the region of these low-latitude stations is approximately 30 000 nT. Figure 2 shows the time-series from each of the stations, with their mean values removed. The most obvious signatures are the diurnal Sq variations that result in minima near local noon; these exist at all stations but are clearest at Hermanus, primarily due to the socalled "coastal effect" which arises due to the magnetic field induced by time varying magnetic fields in the conducting ocean water (Pal'shin et al., 1999).
The D st index over this period indicates slowly varying fields and no magnetic storms although several substorms were observed in the auroral signatures of the CANOPUS photometers and magnetometers. Since D st is computed by the convolution of magnetometer signals from stations at latitudes lower than about 36 • (L < 1.53) it gives a good indication of the state of the inner magnetosphere at the latitudes we are examining. The mean and maximum values of D st during this period investigated were −15 nT and −43 nT, respectively, which represent a relatively small ring current energy density. Table 2. For each station, the columns list the autocorrelation time, τ c ; the power spectral exponent, β; the corresponding Hurst exponent, H β , calculated from β = 2H β + 1; the Hurst exponent, H R/S , calculated by the R/S method; and the fractal dimension, D R/S , calculated from D R/S = 2 − H R/S . The error listed is that from the least-squares regression and so it is only a minimum error as there are other possible sources of error that are difficult to calculate. (The error of D R/S is the same as that for H R/S  have not relied on a single technique but have investigated these properties using two different methods described in the following section, viz. (1) power spectral analysis and (2) range scaling analysis (Hurst, 1951).
A third, more common, method of discerning nonlinear behaviour from a time-series is the use of embedding dimension analysis to evaluate the correlation dimension ν (Grassberger and Procaccia, 1983). The value of ν is determined by counting the number of pairs of points in the time-series that are separated by less than the distance r. This "correlation integral" should scale as r ν for small r. There are two major difficulties with the application of this method. First, the pairs of points must be no nearer than the autocorrelation time τ c (e.g. Hilborn, 1994). For the time-series in this study, τ c ranges from about 4-7 h (these are listed in Table 2) which is a significant fraction of the total length of each timeseries. Even when this restriction is lifted by using only distant points, the strong periodic modulation due to Sq can lead to anomalously large ν (Shan et al., 1991). In the present case, we obtained values for ν ∼ 5 which do not agree with the more robust methods shown below, presumably due to the large correlation times as well as the Sq modulation at 24 h. This is consistent with previous results which indicate that correlation dimensions are adversely affected by strong periodic modulations (Shan et al., 1991), and when the correlation time is large (Shan et al., 1991;Prichard and Price, 1992).

Power spectrum analysis
A great deal of space physics data is self-affine (Ohtani et al., 1995;Takalo et al., 1999) with a power spectral density of the form P (f ) ∝ f −β , where β is the spectral exponent. The power spectrum of the time-series from the Ellis station is shown in Fig. 3 and the best-fit line (calculated over three orders of magnitude in frequency) indicates a spectral exponent of β = 2.093 ± 02.020. All six stations have power spectra  Table 2. With the exception of the Lanseria station, β is slightly larger than 2. This is significant because a value of β = 2 corresponds to a random walk and, as investigated below in Sect. 3.2, a larger value indicates some "persistence" in the time-series (Mandelbrot, 1983).
It appears that β tends to increase with increasing geomagnetic latitude. A linear least-squares fit gives β = (0.95 ± 0.39)L+(0.66±0.63) which indicates a correlation between β and L. The correlation coefficient, however, is r = 0.77 which implies only a weak linear correlation. We are, therefore, unsure of the significance of this result but find that it suggests the possibility that there is some latitudinal dependence in the nonlinear statistics. In addition, the autocorrelation time tends to decrease with increasing geomagnetic latitude. A best-fit line gives τ c = (−8.3 ± 2.1)L + (18.7 ± 3.3) (correlation coefficient r = −0.90). Unfortunately, these six stations have only a narrow spread in L, which means that the strength of these trends is only hinted at with the present data. In fact, we cannot say with certainty that either of these correlations is linear. However, these results do suggest the need for a study over a wider range of low-latitude L-shells which could, quantitatively, discern a trend and perhaps the physical processes of the underlying dynamics with geomagnetic latitude.

Range scaling analysis
Range scaling (R/S) analysis was developed by Hurst (1951) to study time-series whose underlying processes are independent, though not necessarily Gaussian. Here, our timeseries consists of a sequence of measurements of the total magnetic field B(t 0 ), B(t 1 ), . . . B(t M ), where t 0 = 0, t 1 = τ, . . . , t M = Mτ . The time-series is characterized by an exponent, H , which is a quantitative measure of the selfaffinity of the time-series. That is, H relates the typical change in B, B, to the difference in time t by the scaling law B ∼ t H (Mandelbrot, 1983) where H is in the range 0 ≤ H ≤ 1. This is a nonuniform scaling where the shape of the time-series is invariant under a transformation that scales the coordinates differently and is a hallmark of self-affinity. For the usual Brownian motion, which is a stochastic random walk, H = 0.5. Larger values of H indicate some memory or persistence. Smaller values indicate "anti-persistence," which means that the time-series is more volatile and choppy. One method of determining H is to use R/S analysis. For example, R/S analysis was recently applied to the high-latitude AE index and shown to provide a robust estimator of deterministic chaos (Price and Newman, 2001).
Our analysis consisted of taking the raw positive definite time-series x(t) of length M and taking the first differences of the natural logarithm, thus creating a new time-series B(t).
B(t p ) = ln(x(t p+1 )) − ln(x(t p )); p = 1, 2 . . . , M − 1 (1) Following this we take the time-series B(t) and subtract the sample mean B to obtain a new series Next, a cumulative time-series, Y , is derived and an adjusted range, R, is formed in terms of the maximum minus minimum value of the cumulative series Y , R = sup(Y 1 , Y 2 , . . . , Y T ) − inf(Y 1 , Y 2 , . . . , Y T ). The rescaled range, R/S, is then given by the ratio R/σ , where σ is the usual standard deviation. This quantity scales, with respect to T , by the power law where T = t n , and H is the Hurst exponent. The value of H can then be evaluated from a plot of log(R/S) versus log(T ) and a measurement of the slope of the best fit line. Figure 4 shows the rescaled range R/S for the Ellis station, and a best-fit line is shown, resulting in a Hurst exponent of H R/S = 0.504 ± 0.017. The rescaled ranges for the other stations result in similarly good linear fits. The Hurst exponents, H R/S , for all the stations as calculated by the R/S analysis are listed in Table 2. The uncertainties in the calculated values for the Hurst exponent were calculated by starting with the uncertainty in the measured magnetic field, ± 0.1 nT. These errors were propagated through the calculations listed above to obtain the error value for R/S (see Fig. 4). The best-fit slope in Fig. 4 was obtained through linear least-squares analysis taking into account the error in the dependent variable (e.g. Press et al., 1988). All of the linear fits were consistent with the data with chi-square probabilities greater than 0.99, that is, the linear slopes are highly significant. There is no reason to expect, a priori, that there should be a linear relationship between β and H R/S . However, for self-affine data, there is a relation between the Hurst exponent and the spectral exponent, β = 2H + 1 (Turcotte, 1992, p. 78). In Table 2, therefore, the Hurst exponents as calculated from the measured spectral exponent (H β ) are also listed. Although the power spectrum and R/S analysis are independent modes of investigation, these two techniques of calculating the Hurst exponent are relatively consistent although they do not always agree within the estimated errors. We do not expect a perfect correspondence since the methods of investigation are independent and the calculation of R/S is significantly more stable against sudden phase changes of fluctuations than the calculation of a power spectrum. A more rigorous method to measure the strength of their agreement is to determine if the spectral exponent, β, is correlated to the Hurst exponent, H R/S . A linear least-squares fit results in the relation H R/S = (0.54 ± 0.13)β − (0.62 ± 0.29) with a correlation coefficient of r = 0.90. This result leads to the conclusion that these magnetic field time-series are selfaffine (see Fig. 5). In addition, this is evidence that the data behave in a self-organized manner and the calculation of the fractal dimension is meaningful. For a time-series that can be modeled as fractional Brownian motion (Mandelbrot, 1983), the relationship between H and the fractal dimension D is Under these assumptions the fractal dimensions of all stations, as calculated from H R/S and Eq. (5), are listed in Table 2. As expected, the fractal dimensions are all near 1.5 although the stations at the higher latitudes exhibit a somewhat lower value for D; this is consistent with the fact that the Hurst exponents exhibit somewhat more persistence. Similar to the autocorrelation times and spectral exponents, the Hurst exponent exhibits a weak correlation with magnetic latitude or L-shell. The four lower latitude stations are consistent with H ∼ 0.5, indicating stochastic behaviour. On the other hand, the Bos and Herm stations, located at the highest L-values, demonstrate persistent behaviour (as indicated by a value of H > 0.5) which means that these time-series have long memory effects. In the language of nonlinear dynamics, the data exhibit a sensitive dependence on initial conditions, one of the hallmarks of chaos. This feature, the latitudinal dependence of nonlinear features in magnetic time-series, is strong enough to warrant further research. Analysis of data from other magnetometer chains (IMAGE, CANOPUS, MEASURE) is underway.

Possible systematic errors and sources of bias
Because we have used relatively short time-series, it is a reasonable concern that our estimates of the Hurst exponent are affected by the length of the time-series rather than by actual dynamics; a time-series that is too short may bias the estimates. We have investigated this possibility by randomly reorganizing the data so that the order of observations is completely different from that of the original time-series. Because the actual observations remain the same, the frequency distribution of the time-series remains unchanged. If there was a long memory effect in place, the order of the data would be very important and the scrambling effect would be to destroy the structure of the system, thus resulting in a much lower Hurst exponent estimate. However, if the length of the time-series is resulting in bias, scrambling can have the opposite effect, resulting in a Hurst exponent that is even larger than the original estimate (Peters, 1991, p. 75). For stations 5 and 6 which showed persistence, we found that scrambling the original series caused a drop in the value of the Hurst exponents which shows that the long memory process was destroyed by the scrambling process. The other four stations, of course, were already stochastic and the reordering process left their Hurst exponents effectively unchanged.
Another reasonable concern is the affect of the periodic Sq variations which might increase the value of the Hurst exponent in an unphysical manner. We investigated whether this had an effect on our analysis by considering the data from the Hermanus station for the whole of January 1993 (unfortunately, the other stations were only temporarily in operation for the 5 days reported here). Range scaling analysis was performed on a new time-series obtained by subtracting the mean of the three quietest days of the month (21, 22, 23 January) from the original time-series (this essentially removes the Sq effect). There was no significant change in the Hurst exponent. We further investigated this possible effect on an artificial (chaotic) time-series for the Lorenz attractor. The Hurst exponent was calculated for this time-series, and then compared to the Hurst exponent that was computed when the time-series was added to a sinusoidal curve that had 5 periods for the entire length of the series. The two Hurst exponents were statistically equal (i.e. within the error bars).

Discussion and conclusion
Several previous studies have investigated the fractal properties of high-latitude geomagnetic variations through the use of the AE and AL indices. The number of degrees of freedom that describes the system, as measured by the correlation dimension, was found to be between 2.2 and 4.2 (Vassiliadis et al., 1990;Roberts, 1991;Shan et al., 1991;Sharma et al., 1993). In the present paper, low-latitude geomagnetic variations have been examined using different fractal techniques that are not subject to difficulties due to long correlation times and strong periodicities. As stated above, individual magnetometer measurements were used because they have higher spatial and temporal resolution than a global index such as D st .
The L-shells of the stations ranged over L = 1.47 − 1.83, indicating that the stations are sampling a very different region of the magnetosphere, namely the plasmasphere, to the high-latitude stations used to compute AE and AL. The period studied encompassed low dynamic magnetospheric activity, as indicated by relatively steady D st values. This implies that the influence of the ring current perturbations is not large and that the perturbations are primarily due to the ionospheric solar quiet and auroral electrojet currents, as well as possible plasmaspheric influences. The previous studies mentioned have argued that the dynamical processes associated with substorms and, in particular with the AE and AL indices, are low dimensional. Here, we conclude that an even lower dimensional behavior characterizes the low-latitude magnetosphere. The average fractal dimension obtained was near D ∼ 1.5 which is lower than the high-latitude results. This might indicate that the inherent dynamical properties of the low-latitude magnetosphere are less complex and can be described by fewer degrees of freedom than the high-latitude magnetosphere that was previously examined via global indices. Qualitatively, the good agreement between the spectral exponent and the Hurst exponent H R/S lends credence to the existence of self-organized behaviour at low latitudes.
Perhaps the most interesting result is the correlation of every statistical property of the magnetic field time-series with magnetic L-shell. In particular, the increase of H R/S (to values greater than 0.5) with increasing L suggests that the complexity of the inner magnetosphere increases with magnetic latitude and L-shell. However, it is impossible to be certain of the significance of these results since we have only a small sample of data (only six stations over a very narrow latitudinal range). A more detailed global study, with many latitudinally spread stations, is necessary to determine whether the behaviour found here is consistent with that at higher latitudes. Of course these trends, if extrapolated to high latitudes, would lead to unreasonable values in the auroral region. For this reason, we conclude not that the complexity increases linearly with L but only that the trend is suggested and further study is needed.