Science Website Traffic in Earthquakes

This electronic supplement contains tables of Akaike information criterion (AIC) and Omori-law parameters for Heaviside step function and figures of daily seismicity rates, website visits, and comparison of fits of Omori-law, exponential, and power-law distributions to website traffic data.

We obtained earthquake data from GeoNet (see Data and Resources) for the time window reported in Figure 1 in the main article. We used the updated GeoNet earthquake catalog (downloaded in March 2014) for a spatial domain delineated by −43.096° ≤ latitude ≤−43.964° and −171.404° ≤ longitude ≤ −171.907°. On the basis of deviation from *b*-values at *M*_{L} <3.0 in Gutenberg–Richter scaling for the Canterbury earthquake sequence (CES; Shcherbakov *et al.*, 2012), we use *m*_{c} 3.0 in our analysis. In our analysis of seismic data following the 7.1, 6.2, and 6.0 earthquakes, we use 0.1 ≤ *t* ≤ 100 days. We downloaded earthquake magnitude (*M*_{L}) and processed geometric mean peak horizontal ground acceleration (PHA) data for CES earthquakes recorded at the Christchurch Hospital (located in central Christchurch) from GeoNet (see Data and Resources).

We obtained website traffic data for the four websites considered herein from the sources listed in Acknowledgments for the website traffic monitors listed in Table 1 in the main article. The time stamps for website and earthquake data were corrected to UTC to enable derivation of time since mainshock (but the dates referred to in the main article are New Zealand Standard Time [NZST]). Definitions for “Unique Visit(or),” “Total Visits” (or “Visitors” or “Visits”),” “Pages,” “Hits,” “Bandwidth” (or “Data Transfer”) follow the protocols defined at AWStats website (see Data and Resources).

We used the modified Omori law for aftershock decay

$$r(\ge {m}_{\mathrm{c}},t)=dN/dt=K/{(t+p)}^{c},$$

in which the rate of earthquakes (*r*) with a minimum cutoff earthquake Richter magnitude considered in the analysis (*m*_{c}) for a given time interval after the mainshock (*t*) can be calculated using modified Omori law. We used the maximum-likelihood estimation (MLE; Bhattacharya *et al.*, 2011) to derive Omori-law parameters following these events. We obtain ranges of *K* = 41–230, *c* = 0.001–0.1 days, and *p* = 0.8–1.0 (Fig. 3 in the main article).

To make the daily website traffic suitable for use in the MLE, which requires times of individual events (website visits in this case), we generated a simulated distribution of visits per day. For each day, the total number of visits was assigned a time within that day, based on an assumed normal distribution with a mean centered at noon.

To account for pre-earthquake website traffic, we derived mean daily traffic for three months (see Data and Resources) and eight months (see Data and Resources) and subtracted the mean daily pre-Darfield earthquake totals from each of the postearthquake daily totals. To normalize the website data, we divide the daily traffic by the respective background, daily pre-Darfield visit rate for DrQuigs and GeoNet.

To derive a modified version of the Omori law for website traffic, we used the same MLE technique as we did for the earthquakes. To perform the MLE, we used the MLE function within MATLAB to minimize the negative log-likelihood function using the *fmincon* optimization function. We placed no bounds on *K*, other than requiring it to be positive; we conservatively allowed *p* to vary between 0.32 and 1.5; and we allowed *c* to vary between 0 and 30 days. The bounds on *p* and *c* were chosen via iterative running of the optimization routine, progressively relaxing the constraints on *p* and *c* until none of the optimized data resulted in a value of *c* or *p* that were set to the bounds (indicating that the solver was able to find a true minimum of the likelihood function).

Although the use of Omori law to describe aftershock rates is well established, there is limited previous work considering its application to website traffic data, thus it is important to consider whether Omori law is the most appropriate model for our data. As alternatives, we consider two common distributions that share some features with Omori law, an exponential and power law. Although it is common to fit either exponential or power-law functions via procedures like minimization of least-squares error, to directly compare a general exponential or power law to the Omori law we use here, it is necessary to use a similar approach of fitting for all, in this case MLE. To do this, we consider the probability density functions (PDFs) of the general forms of both the exponential and power-law distributions. The PDF of the exponential distribution is defined as

$$p(x)=\frac{1}{\mu}{e}^{\frac{-x}{\mu}},$$

in which *μ* is the mean of the sample data *x*, and the general form of a power-law distribution is simply

$$p(x)\propto {x}^{-\alpha},$$

in which *α* is a scaling parameter, and *x* is the sample data. To estimate the parameters, we use the same input for the sample data as the Omori-law function (see discussion in the previous section). For the exponential distribution, we use the *expfit* function within MATLAB to perform the MLE. For the power-law distribution, we use the methods outlined by Clauset *et al.* (and use the MATLAB code *plfit* provided by the authors, available at Santa Fe Institute website; see Data and Resources). We estimated parameters for the exponential and power-law distributions for both GeoNet and Dr. Quigs web traffic data for the three large events discussed in the main article. Visual comparison of these fits suggests that the Omori law best fits the data (Fig. S1), but we also consider a quantitative estimation of the goodness of fit using the AIC (Akaike, 1974). The AIC is defined as

$$\mathrm{AIC}=2k-2\mathrm{ln}(L),$$

in which *k* is the number of estimated parameters, and *L* is the maximum value of the likelihood function. The AIC can be useful to compare series of models, with the model yielding the lower AIC value being preferred. Calculated AIC values for the three distributions (i.e., Omori law, exponential, and power law) confirm the visual comparison, with the AIC values of the Omori law being significantly lower than AIC values for both the exponential and general power law (see Table S1).

We created the MAG-PHA value presented in the main article by dividing individual *M*_{L} and PHA values by the mean values for the 34 earthquakes with largest PHA values and multiplying the mean-normalized *M*_{L} × mean-normalized PHA to obtain “MAG-PHA.” We calculated the daily change in normalized website visits for 17 of these 34 events spanning 498 days of the earthquake sequence. We use the normalized website traffic because this is more generally applicable; that is, the normalized website traffic can be interpreted as a multiplier of the background web traffic for a given site and is not tied to the particulars of traffic to GeoNet. In several cases, multiple events occurred on the same day, and thus we cannot differentiate how much of a step change in website traffic is due to individual earthquakes on that day because we only have summed daily traffic data. We used the maximum MAG-PHA for a day with multiple events but show the standard deviation of the MAG-PHA for all events on that day in Figure 4a in the main article. For a given event, we calculated the instantaneous change in normalized website traffic by subtracting the normalized website traffic on the day of the event from the traffic on the following day and subtracting the traffic from the previous day from the traffic on the day of the event and then using the maximum of these two quantities. This procedure is to account for the fact that we only have daily traffic, and depending on the timing of the earthquake, the website visit response may happen on the day of the event (e.g., if the earthquake happened early in the morning) or on the following day (e.g., if the earthquake happened late in the day). In the case of two events on consecutive days (typically within hours of each other), as is the case for earthquakes on 21 and 22 February, both events are characterized by the same step change, thus when fitting a relationship between MAG-PHA and change in website traffic, we excluded the 22 February event because this was a clear outlier generated by our calculation method. After removing this outlier, we used a least-squares regression to fit an exponential function to the data. We performed a similar comparison between step change and magnitude (*M*_{L}), this time fitting a linear relation based on the appearance of the data. The same 22 February event was excluded, as were two additional events, the initial Darfield earthquake and an aftershock on 24 December 2010. The Darfield earthquake was omitted because of the novelty effect; the *K*-value was much smaller than predicted by our *K* versus *M*_{L} function because this earthquake was the first in which the utility of GeoNet was just beginning to be realized. Conversely, the 24 December 2010 event generated a much larger *K* for the *M*_{L} because this event was centered right beneath the Christchurch central business district during a slow news day (holiday); we surmise that this resulted in a large *K* response.

Just as with earthquake rates in aftershock sequences, website traffic can be excited by a large aftershock (Fig. 4b in the main article), setting off another Omori-law style decay superimposed on the original. To deconvolve these individual events and to generate a continuous function, we model the website traffic using individual terms for distinct aftershock sequences combined with a heaviside step function. This takes the general form

$$r(t)=\frac{H(t-{t}_{1})\times {K}_{1}}{{({c}_{1}+(t-{t}_{1}))}^{{p}_{1}}}+\frac{H(t-{t}_{n})\times {K}_{n}}{{({c}_{n}+(t-{t}_{n}))}^{{p}_{n}}}+\dots ,$$

in which *r*(*t*) is the normalized daily visits with time, *K _{n}*,

$$H(x)=\{\begin{array}{cc}0,& x<0\\ 1,& x\ge 0\end{array}\mathrm{.}$$

To build this function, we first manually identified the timing of significant events that modify the rate of visits and did not result in single-day spikes in traffic. We used these dates to divide the record of website visits into eight distinct segments. We then fit the first segment (which began with the main event in 3 September 2010 UTC) using MLE. Because the ultimate function is additive across the whole time series, and we are interested in extracting the Omori-law parameters for individual events, we use the parameters from this first fit to detrend the website traffic for all days following the end of the first segment; that is, we calculate what the website traffic would be for the remainder of the record based on this initial decay and subtract those values from the actual values. We then use the MLE technique to find the Omori-law parameters for the second segment on this detrended data and continue this pattern of detrending the following days, using the MLE technique to calculate individual Omori-law parameters (see Table S1 for these parameters and the segment boundaries).

We obtained daily views for a suite of science articles pertaining to different natural disasters from staff at *The Conversation* (see Data and Resources) through an e-mail request. The events and websites are shown in Figure S3. The red line is Omori fit, and the green line is simple power-law fit. Omori parameters and power-law exponent (*α*) are on each graph. Exponential fits were also attempted but did not fit the data well and are not shown in these plots (analysis of exponential fits is available from A. Forte). The general form of the traffic decays adheres to Omori law and a simple power law; in some instances, the former provides a better fit, and in other instances a simple power law provides a better fit.

Table S1. AIC values for Omori-law, exponential, and power-law distributions for the three main events discussed in the main article.

Table S2. Omori-law parameters for Heaviside step function.

Figure S1. (a–c) Plots of daily seismicity rates, (d–f) daily website visits for Dr. Quigs, and (g–i) daily website visits for GeoNet for first 100 days following the Darfield earthquake, 22 February 2011 Christchurch earthquake, and 13 June 2011 Christchurch earthquake. (j)–(l) Website traffic normalized by background rates and earthquake data. Corresponding residuals are shown in (m)–(o). Earthquake and website traffic data are fit by the MLE of Omori law. Values for *K*, *c*, and *p* are shown for each plot, along with percent residuals (see the Methods section in the main article).

Figure S2. Comparison of fits of Omori-law, exponential, and power-law distributions to website traffic data from GeoNet and Dr. Quigs for the 100-day periods following the 4 September 2010, 22 February 2011, and 13 June 2011 earthquakes within the CES. Parameters for fits were estimated using the MLE method for all three models. Fits shown here are identical to those in Figure 3 in the main article.

Figure S3. Comparison of fits of Omori-law and power-law distributions to website traffic data for a variety of natural-disaster-related science stories on *The Conversation* (see Data and Resources). The “simple power-law fit” in the green line is an MLE of a power law, and the red line is an MLE of Omori law. Omori parameters and power-law exponents (*α*) are on each graph.

All of the raw and processed website traffic data presented in this study are freely available upon request (M. C. Q.). Data for websites other than DrQuigs.com may also be obtained by contacting Kevin Fenaughty, Chris Crowe, and Anekant Wandres. The data for this article are also obtained from https://theconversation.com/au, http://www.santafe.edu/~aaronc/powerlaws/, http://www.awstats.org/docs/awstats_glossary.html, and www.geonet.org.nz.

Akaike, H. (1974). A new look at the statistical model identification, *IEEE Trans. Automat. Contr.* **19**, no. 6, 716–723.

Bhattacharya, P., M. Phan, and R. Shcherbakov (2011). Statistical analysis of the 2002 *M*_{w} 7.9 Denali earthquake aftershock sequence, *Bull. Seismol. Soc. Am.* **101**, 2662–2674.

Clauset, A., C. R. Shalizi, and M. E. Newman (2009). Power-law distributions in empirical data, *SIAM Rev.* **51**, no. 4, 661–703.

Shcherbakov, R., M. Nguyen, and M. Quigley (2012). Statistical analysis of the 2010 *M*_{w} 7.1 Darfield earthquake aftershock sequence, *New Zeal. J. Geol. Geophys.* **55**, no. 3, 305–311.

[ Back ]