Changes in the response of the AL Index with solar cycle and epoch within a corotating interaction region

Abstract. We use observations in the solar wind and on the ground to study the interaction of the solar wind and interplanetary magnetic field with Earth's magnetosphere. We find that the type of response depends on the state of the solar wind. Coupling functions change as the properties of the solar wind change. We examine this behavior quantitatively with time dependent linear prediction filters. These filters are determined from ensemble arrays of representative events organized by some characteristic time in the event time series. In our study we have chosen the stream interface at the center of a corotating interaction region as the reference time. To carry out our analysis we have identified 394 stream interfaces in the years 1995–2007. For each interface we have selected ten-day intervals centered on the interface and placed data for the interval in rows of an ensemble array. In this study we use Es the rectified dawn-dusk electric field in gsm coordinates as input and the AL index as output. A selection window of width one day is stepped across the ensemble and for each of the nine available windows all events in a given year (~30) are used to calculate a system impulse response function. A change in the properties of the system as a consequence of changes in the solar wind relative to the reference time will appear as a change in the shape and/or the area of the response function. The analysis shows that typically only 45% of the AL variance is predictable in this manner when filters are constructed from a full year of data. We find that the weakest coupling occurs around the stream interface and the strongest well away from the interface. The interface is the time of peak dynamic pressure and strength of the electric field. We also find that coupling appears to be stronger during recurrent high-speed streams in the declining phase of the solar cycle than it is around solar maximum. These results are consistent with the previous report that both strong driving (Es) and high dynamic pressure (Pdyn) reduce the coupling efficiency. Although the changes appear to be statistically significant their physical cause cannot be uniquely identified because various properties of the solar wind vary systematically through a corotating interaction region. It is also possible that the quality of the propagated solar wind data depends on the state of the solar wind. Finally it is likely that the quality of the AL index during the last solar cycle may affect the results. Despite these limitations our results indicate that the Es-AL coupling function is 50% stronger outside a corotating interaction region than inside.


Introduction
Very early in the space age it was shown that geomagnetic activity is related to the solar wind speed (Snyder et al., 1963) and controlled by the north-south component of the interplanetary magnetic field (Fairfield and Cahill, 1966).This result was interpreted as evidence of magnetic reconnection between the interplanetary magnetic field (IMF) and the earth's dipole field.According to Dungey (1963) the rate per unit length at which southward IMF is transported to the subsolar magnetopause should be proportional to the dawn dusk component of the electric field given by E y =V B z .Subsequent work suggests that the magnetosheath flow pattern and stagnation of the flow may modify this simple assumption.Early work examined the relation of averages of different solar wind parameters versus various magnetic indices finding that Published by Copernicus Publications on behalf of the European Geosciences Union.
the larger is B z , the stronger is magnetic activity (Schatten and Wilcox, 1967).This work also showed that B z in gsm coordinates exhibits the highest correlation with magnetic activity (Hirshberg and Colburn, 1969).Arnoldy (1971) studied the relation between gsm B s (B z northward = 0) and the hourly integral of AE (auroral electrojet index).The hourly integral was calculated as the sum of several samples multiplied by the time between samples (τ ∼approximately 10 m) represented as B s τ .He found that the highest correlation between B s τ and AE occurred when the input was taken one hour ahead of the output, but there was correlation at other lags as well.This led him to express the output AE as a linear combination of input values at lags of 0, 1 and 2 h.The correspondence between the model predictions and the observations was remarkable.This model was actually a linear prediction filter.Meng et al. (1973) used 5-min resolution data to study the cross correlation between AE and IMF B z and found that the peak correlation occurred at ∼40 min delay, a value somewhat less than the value obtained by Arnoldy.The first study to explicitly use linear prediction filters to investigate the relation between IMF B z and various indices was performed by Iyemori et al. (1979).The authors used hourly averages to show that AL, AU, AE and D st are all reasonably well predicted by IMF B z .The auroral electrojet filters were all short with only a few samples contributing to the output, but the D st index depends on inputs for many hours prior to the current time.Clauer et al. (1981) extended this work with 2.5 min resolution data and consideration of three different coupling functions: epsilon, V B s and V 2 B s .Epsilon is proportional to the product of the interplanetary Poynting vector and a "gating function" that depends on the clock angle of the IMF around the Earth-Sun line (Perreault and Akasofu, 1978).All of the coupling functions produced filters that rise rapidly to a peak in an hour or less and then decay more slowly for several hours.They found that the epsilon parameter produced considerably less accurate predictions than the other coupling functions and that its filter was much noisier.They also noted that moderate activity filters tended to peak at about 60 min while strong activity filters peaked near 30 min.They interpreted this as evidence of a possible nonlinear response of the magnetosphere to the solar wind.Clauer et al. (1983) used the same technique to determine the prediction filter relating E s (component of E y due to B s ) to the ring current asymmetry index.They demonstrated that this filter is very similar to the AL filter suggesting a relation between the westward electrojet and the current system responsible for asymmetry in D st .
The nonlinearity of the E s to AL response was investigated in greater detail by Bargatze et al. (1985).The authors selected isolated intervals of activity and then characterized each interval by its median value of AL.Prediction filters were created for successive intervals.They found that the filters consisted of two peaks at about 20 and 60 min.The 60min peak was highest for moderate activity while the 20-min peak was highest in strong activity.This result stimulated a long sequence of papers that utilized this dataset to study solar wind coupling to the westward electrojet.
A number of reviews of solar wind coupling to geomagnetic activity were written about this time.Reiff (1983), Baker (1986), andBaumjohann (1986) described a variety of standard statistical techniques for studying solar wind coupling.Clauer et al. (1986) provided a detailed description of the techniques of linear prediction filtering.McPherron et al. (1988) reviewed results of linear prediction filtering noting that the E s -AL response function can be approximated by a Rayleigh function with time constant one-hour.Since the maximum of a Rayleigh function occurs at a time equal to the time constant this implies that the peak AL response to a delta function input will be delayed by this amount.Note, however, that although this result was quoted in the abstract it was not shown in the body of the published paper.In the results discussed below we demonstrate the truth of this statement.McPherron et al. (1988) also showed that the transfer function is a low pass filter with a cutoff frequency about 0.1 mHz (2.8 h).This filter explained less than 45% of the variance in the dataset suggesting that other factors beside the solar wind are important in the creation of AL.However, the authors showed that more than 90% of the variance in an individual event could be described by a very simple response function consisting of two delta functions of specified amplitude and time delay provided the four parameters are varied from event to event.Their interpretation of this result was that the AL index contains two components: one directly driven by the solar wind through the measured impulse response; and another driven by energy stored in the magnetotail and unloaded in a substorm expansion.The surprising result was that the second component also appears to be directly related to the solar wind electric field.However, the time delay when this component begins is controlled by internal processes and hence on average does not correlate with the solar wind.
Subsequent to the work by McPherron et al. (1988), Weimer (1994) carried out a superposed epoch analysis of AL during the expansion phase and fit a slightly different function to the mean behavior during the expansion and recovery phase.This function was f =c o +c 1 te pt .He found that the decay time constant −1/p decreased from 0.56 to 0.41 h as activity increased.For this function these values imply that the maximum of the response function occurs earlier as activity increases.
Techniques of nonlinear dynamics have been applied to the AE index time series in attempts to identify the type of system represented by solar wind coupling to the auroral electrojets.Vassiliadis et al. (1990) considered the magnetosphere as an autonomous system, one driven by a low-level steady input.In this situation the transient response of the system disappears as the system approaches a semi-steady state governed by internal dynamics.They concluded that the magnetosphere is a low-dimensional chaotic system with a fractal dimension near 4.This implies that only four differential equations are needed to describe the relation of E s to AE.In an extension of this work Vassiliadis et al. (1991) concluded that the Lyapunov exponent of the system, the time to depart from a given state, is only 10 min.Prichard and Price (1992) disputed this result arguing that the conclusion of low-dimensionality and short-lived chaotic behavior was an artifact of a long autocorrelation time in the AE index time series.They concluded instead that the index series usually represents random behavior with some nonlinear structure.Such time series are produced by a driven system with a random forcing function.They suggested that accurate dimensional estimates can only be obtained when the driver (E s in solar wind) is steady for times long compared to the time constants of the system.Price and Prichard (1993) examined one such event and concluded that there was some evidence for deterministic nonlinear coupling.Further analysis of a longer data series by Prichard and Price (1993) again concluded that there is no evidence of low-dimensional chaotic behavior.
Singular spectrum analysis was applied to the AE index time series by Sharma et al. (1993) who again concluded that that the magnetospheric system could be represented by a low dimensional system.Vassiliadis et al. (1993) demonstrated the feasibility of this by representing the system by an LRC circuit.Parameters in this model were determined by least square optimization.When driven by the solar wind these simple low-dimensional models were able to reproduce the behavior in the AE time series about as well as linear prediction filters.Quantitatively the authors showed that averaged over 1-2 d intervals their LRC models usually predicts more than 40% of the AE variance.This is very close to the average value of 45% obtained later in this paper.
The first attempt to use modern techniques of information theory to treat the magnetosphere as an input-output system was carried out by Price et al. (1994).The authors used "local linear filters" to predict the AE index some time ahead of the current time.This technique assumes that the current state of the system is defined by a state vector consisting of a sequence of previous values of the input and output, both normalized by their respective standard deviations.The historical record is searched to find previous examples of this state.An ensemble of these "nearest neighbors" is used to calculate a filter to advance the prediction one time step.Single-step prediction uses the last measured values to advance the state and prediction.Multi-step prediction uses the previously predicted values of the output and measured values of input to advance both the state and prediction.Price et al. (1994) found that their prediction errors do not stabilize until at least 500 nearest neighbors are used to advance the prediction.For single-step prediction about 95% of the variance in the next value is predictable.However, in multi-step prediction the prediction efficiency stabilizes at about 60% after 60 min.Twenty different coupling functions including gsm V B s and epsilon were considered with single-step prediction efficiencies that were virtually identical and close to the result of persistence, i.e. the next values is the same as the last value.The authors conclude that there is little evidence for nonlinear coupling.They support the conclusion of Bargatze et al. (1985) and McPherron et al. (1988) that the AE time series contains a strong and unpredictable stochastic component.
An extension of input-output system analysis has been reported by Vassiliadis et al. (1995).The authors utilize the same basic principles as did Price et al. (1994), but with a number of differences.Among these were the use of the AL index rather than AE which includes AU and AL indices; the use of the solar wind monitor IMP-8 closer to the Earth than ISEE-3; calculation of both moving average and autoregressive moving average filters.In addition they performed detailed optimization of the various parameters used in the analysis method including the length of the filter, the separation of samples used in the filter, the number of nearest neighbors defining the state space, and the number of singular values used in the matrix inversion.The authors find that local moving average (MA) filters that depend only on the current state of the solar wind are optimum for a 2.5 h long filter when the input series is sampled at five minute intervals, with 100 nearest neighbors, and separated from the current state by at least 10 h.Such filters make single step predictions of AL with a prediction efficiency of order 75%.Autoregressive moving average filters (ARMA) that depend on the state of both the solar wind and previous AL index do much better (∼90% of variance) with far fewer coefficients (4-6) in the two parts of the filter.When these filters are iterated for about 4 h using previous predictions of AL and the observed input series they still predict about 65% of the variance.In comparison linear filters predict about 40% of the variance.(Note the authors report prediction quality with correlation coefficients which are approximately the square root of prediction efficiency.)An important point made by the authors is that the AL predictions are stable against perturbations of the initial conditions used to start an iterated prediction.
In a later paper Vassiliadis et al. (1996) these same authors reexamine the question of whether the V B s -AL coupling in nonlinear.They conclude that the answer to this question strongly depends on the details of the analysis.In particular, when filters are averaged over a large range of activity levels they are biased toward becoming linear prediction filters.They conclude that the magnetosphere is nonlinear and that this must be taken into account by the use of state-dependent prediction filters.
Until recently little additional work has been done on the problem of V B s -AL coupling.Attention has turned to other indices such as the PC index and the D st index.Also an effort has been made to define a better input parameter than V B s .For example, Newell et al. (2007) Pulkkinen et al. (2007), however, have investigated the V B s -AL coupling problem using superposed epoch analysis.A set of 150 electrojet activations defined by the onset of a negative bay in AL was selected.The dataset was divided into high and low driving with E y =4 mV/m as the dividing line and high and low dynamic pressure with P dyn =3 nP separating the two classes.They find that the ratio |AL|/E y for low driving before onset is ∼130 and after is ∼180.For high driving the corresponding ratios are ∼110 and ∼140.Similar results were found for high and low dynamic pressure.For low pressure before onset the ratio is ∼130 and after onset it is ∼180.For high pressure the corresponding values are lower, ∼75 and ∼130.Thus weak driving and low dynamic pressure both result in stronger coupling.These results seem to confirm the conclusion of Vassiliadis et al. (1996) that the AL index has a nonlinear response to the solar wind electric field.
The purpose of this paper is to determine whether the nonlinear results obtained by Pulkkinen et al. (2007) are evident in the properties of linear prediction filters.For this analysis we use a technique somewhat like the first method described by Vassiliadis et al. (1995) (see preceding discussion).These authors used a state vector depending only on the solar wind input (E s ) to define the system state and for each state constructed a moving average filter from a number of very similar states.In our case we define the state of the solar wind differently.In particular we use the stream interface in a corotating interaction region (CIR) as a reference time.We assume that all CIRs exhibit similar values of solar wind parameters as a function of time relative to this reference time.Thus we expect that the linear prediction filter that transforms E s into AL is the same for all CIRs but varies with epoch time.We will show that there is a significant change with epoch time with the weakest coupling occurring at the stream interface when the solar wind electric field and dynamic pressure are strongest.
To illustrate why we might expect a change in the V B s -AL coupling function during the passage of a CIR we briefly review the characteristics of corotating interaction regions.A CIR is formed when slow speed solar wind from one longitude on the Sun is overtaken by high-speed wind from a following longitude.The high-speed plasma can not penetrate the slow speed plasma because of the imbedded magnetic fields.Consequently it compresses the plasma and magnetic field near the interface.This creates a spiral shaped interaction region between the two solar wind streams with the interface between the two streams at the center.Total pressure in the plasma reaches a peak along the interface with a gradient away from ridge of high pressure.Ahead of the interface the pressure gradient deflects the solar wind to the west of the Earth-Sun line and behind it deflects the plasma toward the east.With time the region of elevated pressure propagates away from the interface broadening the interaction region.Inside the leading edge of the CIR the slow wind is accelerated while behind the fast wind is decelerated.In a frame of reference moving with the interface the solar wind flow on the two sides is tangential to the interface.The high-speed stream behind the interface contains large amplitude Alfvén waves propagating outward from the Sun.These waves rotate the IMF southward antiparallel to the Earth's magnetic field enabling dayside magnetic reconnection.The reconnection drives magnetospheric convection which in turn drives field-aligned current closing through the ionosphere.It is the Hall current produced by this closure that is evident in the AE indices.Since the electric field in the solar wind is the rate of flux transport per unit length high speed wind creates a stronger electric field than low-speed wind.Combined with the fluctuations in B z caused by the Alfvén waves it is expected that geomagnetic activity is stronger after the interface than it is before.
An illustration of the average properties of a CIR derived from superposed epoch analysis is presented in Fig. 1.Zero time in this analysis is the time a stream interface, defined as the zero crossing of the azimuthal flow angle of the solar wind passed the Earth.The right side of the figure displays important parameters derived from solar wind plasma and magnetic field measurements.In each panel the shaded regions bounded by the upper and lower quartiles define the range within which 50% of the data fall.The top panel shows that solar wind dynamic pressure begins to increase one day before the stream interface.This is the leading edge of the compression region.Dynamic pressure peaks at the stream interface then decays over the following two days.The behavior of total pressure in the solar wind frame (panel 2) is very similar although starting its increase a little later.The solar wind electric field (E y ) (panel 3) begins to increase 12 h before the interface, peaks at the interface, and decays slowly over a period of three days.Beta (panel 4), the ratio of thermal to magnetic pressure in the solar wind is nearly 2.0 just before the interface, falls to about 1.0 at the interface, and takes several days to return to normal.Solar wind Mach number (panel 5) is typically around 8.0 but inside the CIR it falls to 6.0 and then recovers over the next three days.
The left side of Fig. 1 shows measures of geomagnetic activity during the passage of a CIR.All panels show that activity indices begin to increase a few hours before the CIR stream interface, peak 6-12 h after the interface, and then take about five days to decay to the quiet levels present before the interface.The analysis presented in this paper utilizes E s , the rectified version of E y , as input, and the AL index (a component of AE) as output.If coupling is stronger for low dynamic pressure and a weak driver we would expect to find that the ratio between AL and E s is largest at the edges of the figure and weakest at the center.We will demonstrate that this is the case.

Datasets and preprocessing
The output parameter used in this study is the lower auroral electrojet index (AL) calculated and distributed by the World Data Center -C2 in Kyoto, Japan.We have downloaded these data and interactively edited the AL time series to flag obvious spikes in the index.The input parameter is the rectified dawn-dusk electric field of the solar wind in gsm coordinates.To obtain this quantity we have processed all available solar wind and IMF data from the Wind and ACE spacecraft.The intervals covered by these spacecraft are Wind 1995 to present and ACE from 1 February 1998 to present.
Rectified electric field is calculated from solar wind speed and gsm B z propagated to the subsolar bow shock by the Modified Minimum Variance Method (Bargatze et al., 2005;Weimer et al., 2003).This method uses a moving window to calculate the time-varying normal to discontinuities in the interplanetary magnetic field (IMF).For each window the mean field is calculated and a minimum variance analysis performed on the projection of the field perpendicular to the mean field.The minimum variance direction is taken as the normal to discontinuities in the window.Time delays are calculated by projecting the satellite position vector and solar wind velocity onto this normal.The calculated time delays are added to the times of sequential points.When fast solar wind follows slow wind some parcels of plasma will appear to overtake and pass slower parcels ahead of them.Of course this can not happen in the solar wind because the magnetic field is frozen in the plasma.Instead the plasma and magnetic field near the gradient in speed is compressed.In the usual propagation technique this situation is handled by simply eliminating the overtaken parcels.When the wind speed is decreasing fast particles leave slower ones behind.At the subsolar bow shock the original equally spaced time series is distorted into a time sequence with variable time delay between samples.This sequence is interpolated to the original grid of one-minute samples.If the normal is poorly determined (eigenvalues are nearly equal) or the normal is close to orthogonal to the Earth-Sun line the normal can not be defined or is meaningless.In these cases the time delay is interpolated from adjacent values.
A combined one-minute solar wind and IMF dataset with properties somewhat similar to ours is available from the NASA National Space Science Data Center at the url http: //cdaweb.gsfc.nasa.gov/istppublic/.
Unfortunately we were not able to use this data set for our analysis because of the high density of flags in the output.These data have been propagated with very conservative constraints on acceptable output.In addition to the eigenvalues and normal direction flags the data are flagged whenever any parcel overtakes another parcel.Because of these constraints it is rare to have a flag free interval longer than the duration of an E s to AL filter.We note, however, this characteristic of the NSSDC dataset is not a problem in the generation of the dynamic cumulative distribution functions discussed below.

Linear prediction filters
In this work we utilize linear prediction filters (Weiner, 1942) to study the relation between the solar wind electric field and the AL index.Prediction filters represent the most general linear relation between two time series.In an ordinary linear regression the output of a system at a given time is represented as the sum of a constant and a fixed multiple of the input at that time.A finite impulse filter differs only in the assumption that the output at one time is given by the sum of multiples of the input at that time and earlier times.Since a specific previous output depends on previous inputs it is also possible to represent the system by a sum of multiples of previous output and previous inputs.This latter representation is generally more compact having fewer multiplicative coefficients.In the first case the filter is called either a moving average (MA) or finite impulse response (FIR) filter.In the latter case it is called an autoregressive moving average (ARMA) or infinite impulse response filter (IIR).For simplicity we use moving average filters in this work.
The phrase "impulse response" means that a plot of the filter versus time lag is the output that would be generated when the system is stimulated by a single pulse of unit amplitude.For example, we will show that the impulse response relating the rectified solar wind electric field to the AL index is a miniature negative bay with peak amplitude ∼1.5 nT/(mV/m) and duration of ∼2.5 h.(In our figures we have inverted this response to have positive area.)The Fourier transform of the impulse response is the system transfer function or frequency response.For the lower auroral electrojet index (AL) the transfer function is a low pass filter.
The mathematical representation of an ARMA filter is shown in Eq. (1).
To obtain a moving average filter the autoregressive coefficients a i are all set to zero except a 1 which is set to 1.0.In this case the output at the n-th sample point (time) is a sum of multiples (b i ) of inputs I n−i+1 at earlier time lags (i).(This is a convolution.)The set of equations obtained by allowing the index n to take on many successive values can be represented as illustrated in Eq. ( 2).
The left hand side of this set of equations is a column vector of length N+1 corresponding to the segment of the output time series beginning with sample n and ending with sample n+N.The right hand side is the result of a matrix multiplication between a rectangular matrix (X) with N+1 rows and N b columns and a column vector of N b unknown filter coefficients.This set of equations has a simple matrix representation shown in Eq. ( 3).
The first column of the matrix (X) is the input time series corresponding to the output series on the left hand side.The second column is the same series shifted down by one sample, the third column is the input shifted down two samples, and so on until the last column which is the input shifted down by N b +1 samples.This matrix is often called the "design matrix".To determine all N b coefficients the design matrix must have at least N b rows.Usually there are many more rows than coefficients (N >N b ) and the coefficients are over determined by the data.The least square solution for the coefficient vector b is obtained by multiplying both sides by the transpose of the design matrix and then inverting the resulting square matrix as shown in Eq. ( 4).
It can be shown that the product matrix (X T X) has rows and columns that represent the autocorrelation function of the input at various lags.Similarly the product X T O is the cross correlation between input and output as a function of lag.
Often this solution vector is too noisy to be useful and a different solution method is required.An alternative solution technique is called singular value decomposition or SVD (Press et al., 1986).It is shown in matrix algebra than any real rectangular matrix can be represented as the product of three special matrices X= (U) (S) V T .Matrix U is column orthogonal, matrix S is diagonal sorted descending, and matrix V is fully orthogonal.We use these facts to solve for the coefficient vector obtaining Eq. ( 5).
The reciprocal of the diagonal matrix (1/S) is also diagonal with elements in ascending order.Very small elements along the original diagonal become very large elements in the reciprocal matrix.These large elements are the source of noise in the least square solution.The secret of SVD is to set all elements beyond a certain singular value to zero in the reciprocal matrix.This eliminates the terms causing noise in the solution.
This procedure was modified to allow for the possibility of acausal filters.An acausal filter is one in which there is an output before the input is applied.This situation arises in our calculations if the solar wind has not been properly propagated from the upstream measurement point to the magnetopause.The required modification is to time shift the input data with both positive and negative lags.In our case we typically used 60 min before the expected output and 180 min afterwards for a total of 240 lags.

Ensemble matrices
To determine the impulse response for a particular state of the solar wind we must select a number of examples of the particular state and create an average filter for that state.We then systematically vary some parameter such as solar wind dynamic pressure and for each new range of values calculate a new filter from an ensemble of events satisfying these conditions.In this work we have used the Earth passage of a corotating interaction region (CIR) to establish a sequence of states and have calculated the filter relating the solar wind electric field to AL as a function of time relative to the stream interface at the center of the CIR.We then examine this sequence of filters and determine whether there is an observable change in the properties of the filter with epoch time.
The underlying assumption is that all CIRs create similar solar wind states at the same location in a CIR.
We begin by identifying all interfaces between low-speed and high-speed solar wind in the years 1995-2006(McPherron, et al., 2008a, b), b).For each interface we selected a 10-day interval centered on the interface and extracted a segment of solar wind or index data from our original database.The intervals were placed in the rows of an ensemble matrix.Spline interpolation was used to eliminate short gaps (ten minutes or less) in the original data.For the V B s -to-AL filter we constructed matrices of solar wind speed, interplanetary magnetic field (IMF) B z in gsm coordinates, and AL index.The matrices for solar wind speed (V) and B z were multiplied element by element to construct a new ensemble of V B s (rectified (V B z )).The matrices V B s and AL were then used as the input and output of the magnetospheric system.Note that southward B z (B z <0) produces negative E s (as defined here) so that the impulse response between E s and AL is a mostly positive curve.
To calculate the ensemble average impulse response we utilize 24-h sequences of data centered at 00:00 UT on each day around the stream interface.The design matrix for a given day of data and specific stream interface was constructed from the appropriate row in the input ensemble array.For day −5 (left end of a row) the time shifting required to construct the design matrix introduces flags (missing data) from outside the interval.The remaining days simply shift data from the preceding day into the matrix.For acausal filters flags are also shifted into the matrix from the right end of the arrays.In some cases, even after interpolation, there are missing data flags at arbitrary locations in the original data.These flags are shifted downward introducing flags in succeeding rows of the design matrix for that day.At the same time we construct the vector of system output (AL).As a final step we horizontally concatenate the output vector and the design matrix and search each row of this concatenation for missing data flags.Any row containing one or more flags is eliminated from the final output vector and design matrix.In a few cases there were too few rows remaining to calculate a prediction filter.
At this point the ensemble average prediction filter was calculated by two different methods.In the first we used the data available for a given day relative to a stream interface and SVD analysis to calculate a daily filter.These filters were then averaged over all interfaces in a year producing an ensemble average filter for the year.Filters calculated in this manner were highly variable from day to day, but generally represent the data from which they are calculated with high prediction efficiency.Prediction efficiency is defined as shown in Eq. ( 6).The operator "Var" means variance, i.e. mean square deviation of the argument from its mean.
In the second method used for this study we vertically concatenate the design matrix for each day into one long matrix for an entire year.SVD was then used to invert this matrix.Typically this matrix had about 365 000 rows and 240 columns.The filter obtained in this manner is the optimum linear representation of the relation between input and output data for CIR events during a given year.Since there are many different events in a year the prediction efficiency determined for an entire year is lower than the average prediction efficiency of the daily filter calculated in the first method.Typical values of prediction efficiency are about 45% of the variance with an annual average filter.Single day filters typically represent 65% of the variance.

Characterization of filters
In an ordinary linear regression a single number that characterizes the relation between the input and output is the multiplicative constant.In a linear filter the equivalent quantity is the area under the impulse response function.This can be clearly seen by considering an example of a constant input of unit amplitude to the system.After the transient response time of the filter (the length of the filter) the output stabilizes at a value given by the sum of products of the magnitude of the input filter and the coefficients.Because the input is constant it can be taken outside the sum in Eq. ( 1) leaving the sum of the filter coefficients.In this work we use the area under the response function to quantify the response of the magnetosphere.If the coefficients change in a systematic way with the state of the system we conclude that the system is nonlinear.For stream interfaces we use time relative to a stream interface (epoch time) as the state variable.

Example data and filter
An illustration of the relation between rectified solar wind dawn-dusk electric field and the AL index is presented in Fig. 2 for several hours in January 1999.The electric field in the top red trace is zero when IMF B z is positive.During these times there is no variation in the output.Note that in this analysis we removed mean values for a month prior to calculation of the prediction filters and then add the means back so slight differences in the baseline of the observed and predicted values may be present on a given day.When B z is negative E s =−V B s takes on positive values and AL responds after a short time.The predicted response follows the general pattern of the observations but is much smoother, does not contain the extreme variations in AL, and often shows timing differences relative to the observations.These differences cause the prediction efficiency to be significantly less than 1.0.The heavy dashed line (black) in the upper panel is the result of fitting an offset Rayleigh function to the impulse response.This function is given by the equation In this example the three model parameters have the values [−0.0699, 0.0629, 41.9291].The amplitude constants have units of nT/(mV/m) while the time constant is in minutes.For a Rayleigh function the peak response occurs at a time equal to the time constant.After the peak the response dies slowly passing through zero at a lag of 120 min.The response beyond this lag is not usually present in monthly filters and is probably an artifact of detrending by only subtracting the mean.This filter was characterized by the area under the filter between zero and 180 min.We also tried to fit this response function with the Weimer (1994) function discussed in the introduction.We found that this function fits the slope of the rise of the filter somewhat better, but it overshoots the peak value by almost a factor of 2. During the decay this model compensates for the overshoot by being less than required to fit the prediction filter.
The transfer function (frequency response) is defined as the Fourier transform of the impulse response.We have calculated this for each of the filters shown in the top panel except we have not plotted the function for the least square solution which is quite noisy.The heavy dashed black line is the transform of the fitted Rayleigh function.The heavy red line is the solution retaining the largest number of singular values (30).Other traces shown by thin blue lines correspond to the use of fewer singular values in the matrix inversion.It is evident how the use of singular value decomposition suppresses noise (and high frequency components) in the impulse response function.
The transfer function for the V B s -AL impulse response function is a low pass filter.Two arrows in the figure represent important cutoffs.The first at low frequency is physically meaningful and represents the behavior of the magnetosphere.Electric field fluctuations with periods shorter than ∼4 h are attenuated while those with longer periods are not.The second cutoff at higher frequency depends on the number of singular values retained in the solution.The cutoff at ∼0.4 h corresponds to 30 singular values.Solutions using fewer singular values have lower frequency cutoffs.
We also examined the phase response of this filter (not plotted).From zero to 0.2 mHz the response is linear and given by the function φ=af where a=−1.6148×10 4 radians/Hz.This translates to a uniform time delay of 42 min for these low frequency signals.

Properties of interplanetary electric field and AL index relative to CIR stream interfaces
Our emphasis in this work is an investigation of whether prediction filters for magnetic indices change systematically relative to a CIR stream interface.In work reported elsewhere (McPherron et al., 2008a) we have determined 394 stream interfaces in the interval 1995-2007.The behavior of the dawn-dusk component of the gsm interplanetary electric field E y =V B z is shown in Fig. 4.During quiet times several days before the interface roughly half of all electric field values have magnitude less than 1 mV/m.About 12 h before the interface the electric field strength begins to increase.It reaches a value twice as high at the interface.Subsequently it decays slowly reaching typical background strength after three days.
The AL index mimics its driver E y as can be seen in Fig. 5. AL begins to decrease 12-24 h before the interface and reaches a minimum value about 8 h after.The minimum median value of AL is only −150 nT.It is more negative than −300 nT less than 25% of the time.Note that activity is elevated for more than five days after the interface even though E y reaches background values in only three days.

E s -AL filters as function of CIR epoch
Prediction filters relating rectified solar wind electric field to the AL index are plotted in Fig. 6.Each trace in this figure represents an ensemble average for a particular 24-h period around the CIR stream interface.Each filter is computed as an average of that 24-h period preceding or succeeding the zero epoch time at the stream interface averaged over events observed during years 1995-2007.However, the 13-year interval does not include 1996 because no AE indices are available this year.Ensemble average filters were determined first for a given 24-h epoch and all data for events observed during   a given year, and then filters from successive years were averaged.All filters have essentially the same shape starting from zero at zero lag, rising to a peak at ∼20 min, decaying slowly to zero at 150 min, and then slightly overshooting with negative values.The filters are virtually identical except for the filter highlighted with red error bars, which shows the filter for the zero epoch time during the day of the stream interface crossing.The error bars on the interface filter are the  The top panel shows that the 3-day interval beginning two days before and including the interface exhibits the highest predictability of about 45%.The day after the interface during the high-speed stream this value drops to its lowest value of 35% and then rises slowly to a value of 42%, somewhat less than the values before the interface.The bottom panel shows that the area under the interface filter ( 102) is significantly less than the values in the days before and after (125 and 112).The values of the filter area adjacent to the interface are slightly lower than those two days and more before or after.Even though the data are more predictable on the day of the interface than any other day, the area under the filter is smaller.This suggests that the solar wind coupling becomes less efficient as a CIR passes the Earth.Note that the bottom panel indicates that the filters after the interface have a more variable area that they do before, a result consistent with the decreased prediction efficiency after the interface.

E s -AL filters as function of solar cycle
Since we have determined filters for each epoch of the CIR and for each year in the solar cycle we can also display average filters versus phase of the solar cycle.Figure 8 displays prediction filters averaged over all CIR epoch times in each year.Each trace in the figure shows the average filter for a given year.Two traces have been highlighted and annotated.The red trace with the largest area ( 163) is for the year 1995 late in the declining phase of solar cycle #22.The black trace with the smallest area ( 88) is for the year 2001 close to the maximum of the solar cycle.
The variations of the filter properties, prediction efficiency and area, with solar cycle are presented in Fig. 9.As shown in the top panel the annual filters typically predict about 45% of the variance, but in 2001 this dropped to a low of 35%.The data suggest a solar cycle variation in prediction efficiency with the highest efficiency in the declining phase of the solar cycle and the lowest near solar maximum.The bottom panel shows the area of the impulse response function versus phase of the solar cycle.The traces suggest a solar cycle effect with the strongest coupling in the declining phase (1995 and 2007) and the lowest around solar maximum in 2001.

Discussion and conclusions
Previous studies of the relation between the solar wind electric field and AL index have used linear prediction filters, local linear filters, and neural networks.In general the nonlinear models obtain the highest prediction efficiency, and  the linear prediction filters the lowest.Clearly the relation between the two quantities is not linear.Two types of local linear filters have been utilized to approximate the nonlinear behavior.The most general defines the state of the magnetospheric system by a state vector constructed from the input and output time series.Filters are created from an ensemble of similar previous states at each time in the series being represented by the filter.In this case the prediction filter is continuously varying and adapts to represent the nonlinear system.A second type of filter uses the data immediately prior to the prediction point to define a filter to advance the prediction.Our work uses a variant of these methods.We assume that the solar wind establishes the state of the solar wind and that a linear filter can represent the input-output relation for each state.In our case we assumed that a corotating interaction region (CIR) establishes different magnetospheric states as a function of time relative to the stream interface.We then determined an average filter for each day for a 10-day period centered at the stream interface.The filters were determined using an ensemble of CIRs recorded during a given year.The collection of filters allows us to investigate whether there is solar cycle variation in solar wind coupling or a variation with time in a CIR.
Our results show that there are significant differences in the E s -AL coupling function both with epoch within the CIR and with solar cycle.First the prediction efficiency is reduced in the high-speed stream after the interface as compared with the value before and at the interface.Second the area under the impulse response is lowest at the stream interface, slightly stronger the day before and after, and strongest on other epoch days.We find that there is little change in prediction efficiency with the solar cycle averaging around 45% The physical interpretation of these results is not obvious.As shown in Fig. 1 we know that the state of the solar wind changes during a CIR.Two days before a stream interface the density begins to increase as the solar wind speed reaches a minimum value.Twelve hours before the interface the solar wind speed begins to increase reaching a maximum about 1.5 days after the interface.Thereafter the speed falls slowly.The IMF magnitude begins to increase one day before the interface, peaks at the interface, and then decreases over a period of about three days.Temperature of the solar wind is very low 12 h before the interface, a maximum at the interface, and then decays slowly for at least five days after the interface.
These properties of the solar wind affect various derived quantities important to magnetic reconnection.The dawndusk electric field is very asymmetric about the interface as illustrated in Fig. 4. The combination of increased speed and strong magnetic field cause E y to maximize at the interface and to decay slowly as the field strength and speed decrease in the high-speed stream.The density increase before the interface and the speed increase after cause the dynamic pressure to peak at the interface.Plasma beta decreases rapidly from 2.5 beginning 12 h before the interface to 1.5 for many days after the interface.Alfvén Mach number of the solar wind behaves in a somewhat similar manner dropping from a value of 10 to about 7 at the interface.It recovers more quickly returning to normal values in about two days.
Our analysis procedure does not allow us to separate various possible causes of a change in coupling.Is the observed decrease in coupling around the CIR stream interface a result of a smaller magnetosphere due to the enhanced dynamic pressure?Alternatively does the reduced Alfvén Mach number in the solar wind alter the efficiency of reconnection?The superposed epoch analysis reported by Pulkkinen et al. (2007) indicates that either high dynamic pressure or strong driving by E y can reduce the coupling efficiency.To answer these questions we would need to use a different protocol.In particular we would need to create an ensemble of events in which each row is characterized by reasonably constant values of some parameter.We could then calculate impulse response functions for rows selected for a given range of the possible control parameter.Conceivably we could bin according to two different parameters and still have enough data to define the impulse response.To do this we must create ensembles with rows somewhat longer than the response functions, i.e. longer than four hours.This procedure will be difficult and time consuming.The fact that we find significant differences in coupling near the stream interface gives us confidence that the binning procedure would produce interesting results.
We also found apparent variations in E y -AL coupling with solar cycle.The strongest coupling seems to occur in the declining phase in association with recurrent highspeed streams.The weakest coupling appears to be around solar maximum.One possible explanation is that coronal mass ejections (CME) at solar maximum contain solar wind plasma with properties quite different than found in CIRs and these properties affect the size of the magnetosphere and the reconnection process.Magnetic field strength in CMEs is often higher, density and temperatures lower than in CIRs.Another possibility is that the solar wind propagation algorithm works better during the moderate conditions associated with CIRs and the decreased filter area and prediction efficiency near solar maximum is caused by a poor representation of the solar wind arriving at the magnetopause.A third possibility is that the AL index is of lower quality around the last solar maximum (1997)(1998)(1999)(2000)(2001)(2002).This is a time when magnetometer data acquisition from the Siberian sector was very poor.Stations missing from the AE network will reduce the area of the impulse response.Data spikes associated with poor transmission and recording will decrease the prediction efficiency.
Unfortunately, no continuous high quality solar wind data are available before 1995 and insufficient good AE index data are available until the last few years.Thus we can not yet reach a definitive conclusion on whether the solar cycle changes we observe are real effects of the solar cycle, or artifacts of the propagation algorithm, or artifacts of the index generation.In future work we intend to apply this procedure to the Thule Polar Cap (PC) index which is based on one station and has been continuously calculated since 1975 with high quality data and a consistent procedure.If similar trends are apparent in these data we will be more confident that our AL results are physically meaningful.
The application of prediction filters or neural networks to the study of solar wind coupling is based on the assumption that there is a deterministic relation between the input and output variables, in our case V B s and AL index.There is no reason to expect that all of the AL variance is directly controlled by V B s .Some of the variance is likely to be caused by the viscous interaction which presumably does not depend on electric field.Some may be caused by changes in dynamic pressure.Also internal dynamics of the magnetosphere such as substorms in the tail and electron precipitation into the auroral oval will cause changes in the westward electrojet that may have only a probabilistic relation to the solar wind electric field.As explained in the description of the analysis procedure, linear filters depend on the existence of a fixed correlation between the input and output variables.Only that portion of AL that is correlated with V B s is predictable with our analysis.It is likely that multi channel prediction filters would predict somewhat more of the AL variance than a single channel filter.The main result obtained in this work is that the coupling function relating V B s and AL changes significantly with the state of the solar wind as defined by time relative to the interface between low-speed and high-speed solar wind streams.This coupling is weakest when the dynamic pressure and electric field are strongest inside the CIR.A secondary result suggests there is also a solar cycle variation in coupling that is consistent with the result obtained for CIRs, i.e. the coupling is weakest near solar maximum when both dynamic pressure and electric field reach larger values than they do in CIRs.These results show that the V B s -AL coupling function is nonlinear as it changes significantly with the state of the solar wind.Further studies of the role of other variables such as dynamic pressure will provide additional insight into the mechanisms of solar wind coupling to the westward electrojet.

Figure 1 .
Figure1.Results of superposed epoch analysis of solar wind and magnetospheric parameters are plotted versus time relative to a stream interface at the center of a corotating interaction region (CIR).The five panels on the right respectively show solar wind dynamic pressure, total of thermal and magnetic pressure, gsm Ey, plasma beta, and plasma Mach number.The five panels on the left show the AE index, the Sym-H index, the Asym-H index, the 3-h ap index and the PC index.Heavy black lines define the range within which 50% of the data lie.The heavy red line defines the parameter median at each epoch.

Fig. 1 .
Fig.1.Results of superposed epoch analysis of solar wind and magnetospheric parameters are plotted versus time relative to a stream interface at the center of a corotating interaction region (CIR).The five panels on the right respectively show solar wind dynamic pressure, total of thermal and magnetic pressure, gsm E y , plasma beta, and plasma Mach number.The five panels on the left show the AE index, the Sym-H index, the Asym-H index, the 3-h A p index and the PC index.Heavy black lines at the edges of the shaded regions define the range within which 50% of the data lie.The heavy red line at the center of shaded regions defines the parameter median at each epoch.

Figure 2 .Fig. 2 .
Figure 2.An illustration of the relation between time series of solar wind electric field (VBs) and the AL index for two days in January 1999 is presented.Top trace (red) is the rectified electric field.The bottom trace denoted by a thin blue line is the observed AL index for the same two days.The thicker black line is the output of a linear filter representing the relation between these quantities.

Figure 3 .
Figure 3.The VBs-AL prediction filter for the month of January and its frequency response are plotted in the upper and lower panels respectively.The heavy dashed lines show the fit of a Rayleigh function to the filter and its transform (see text for details).

Fig. 3 .
Fig.3.The V B s -AL prediction filter for the month of January 1999 and its frequency response are plotted in the upper and lower panels respectively.The heavy dashed lines in the two panels show the fit of a Rayleigh function to the filter and its transform (see text for details).

Figure 4 .Fig. 4 .
Figure 4.A dynamic cumulative distribution function (cdf) for the gsm IMF dawn-dusk electric field (VBz) for nine years of data measured at the ACE spacecraft.Heavy lines are quartiles of the distribution and thin lines are deciles.The vertical dashed line at zero epoch time is the time a stream interface at the center of a corotating interaction region (CIR) passed the spacecraft.We have reversed the sign of Ey to correspond to Bz.

Figure 5 .
Figure 5.A dynamic cdf of the AL index relative to stream interface (vertical dashed line) within a CIR is presented.Heavy lines are the quartiles as a function of epoch time.Note that the minimum AL index is more negative than -300 nT about 25 % of the time.

Fig. 5 .Figure 6 .
Fig. 5.A dynamic cdf of the AL index relative to stream interface (vertical dashed line) within a CIR is presented.Heavy lines are the quartiles as a function of epoch time.Note that the minimum AL index is more negative than −300 nT about 25% of the time.

Fig. 6 .
Fig. 6.Ensemble average linear prediction filters relating V B s to AL index during the interval 1995-2007.Each trace shows a different day relative to a stream interface.The heavy red line shows the error of the mean of all filters on the filter spanning the CIR stream interface.The filter has been inverted for display purposes.

Figure 7 .
Figure 7. Results for13 annual ensemble average filters are plotted versus time relative to CIR stream interfaces in the years 1995-2007.Blue lines show the annual averages at each daily epoch while thick red lines are the means of all annual averages.Error bars are the standard error of the mean at each time.Note that the day-long intervals used in analysis were centered on the beginning of days.

Fig. 7 .
Fig. 7. Results for 13 annual ensemble average filters are plotted versus time relative to CIR stream interfaces in the years 1995-2007.Blue lines show the annual averages at each daily epoch while thick red lines are the means of all annual averages.Error bars are the standard error of the mean at each time.Note that the day-long intervals used in analysis were centered on the beginning of days.
Figure 8. Ensemble average linear prediction filters relating VBs to AL index during the interval 1995-2007.Each trace is an average over all CIR epochs for a full year.The heavy red line at the top is for the declining phase of the solar cycle for the year 1995.The heavy black line at the bottom is from solar maximum in 2001.

Fig. 8 .
Fig. 8. Ensemble average linear prediction filters relating V B s to AL index during the interval 1995-2007.Each trace is an average over all CIR epochs for a full year.The heavy red line at the top is for the declining phase of the solar cycle for the year 1995.The heavy black line at the bottom is from solar maximum in 2001.

Figure 9 .
Figure 9.The prediction efficiency and filter area of annual average filters for each daylong epoch relative to CIR stream interface are shown by blue circles.The heavy red lines with error bars show the average over all epochs and the standard error of the mean.No AE index data are available for 1996.

Fig. 9 .
Fig. 9.The prediction efficiency and filter area of annual average filters for each day-long epoch relative to CIR stream interface are shown by blue circles.The heavy red lines with error bars show the average over all epochs and the standard error of the mean.No AE index data are available for 1996.
www.ann-geophys.net/27/3165/2009/Ann.Geophys., 27, 3165-3178, 2009 except for the years 2001-2002.The area under the response function is lowest during the rising phase of the solar cycle and again just after maximum.The largest areas are found near the end of the declining phases in 1995 and 2007.