## 3 Method: Empirical Approach  
The most common approach to describing the earthquake occurrence probability is through the survival probability function, that is, the probability that an inter-event time interval will be greater than a given value \(\Delta t_{i}\). This is tantamount to evaluating \(SP(\Delta t>\Delta t_{i})\). In the case of a Poissonian distribution, we have  
\[SP(\Delta t>\Delta t_{i})=exp(-\lambda\Delta t_{i}), \tag{1}\]  
where \(1/\lambda\) is the mean. This means that we can normalize, by simply dividing each \(\Delta t_{i}\) by \(\overline{\Delta t}\). When we have a set of real data \(\Delta t_{i}\) with \(i=1,n\), the observed \(SP\) is simply  
\[SP(\Delta t_{i})=(n-i)/n, \tag{2}\]  
which can also be normalized. The first to propose the usage of real data, instead of trying to model them, was Reference [3], as stated above. I simply extend his approach. There is a further step I introduce before looking at the real data. Even if we do not know what the survival probability function is--and we will likely never be able to model it simply--we can evaluate how far it departs from a Poissonian distribution through a transformation. This means that we can express it as  
\[SP(\Delta t_{i})=exp(-\lambda(\Delta t_{i})\Delta t_{i}), \tag{3}\]  
and, then, define \(\lambda(\Delta t)\) as  
\[\lambda(\Delta t_{i})=-log[SP(\Delta t_{i})]/\Delta t_{i}. \tag{4}\]  
We may immediately evaluate this function if it can be approximated by a Poissonian decay in any \(\Delta t\) interval, as, in that case, it will yield a flattish trend. I believe, moreover, that this approach is the optimal way of evaluating a survival probability function, as we know that it will decrease from _one_ to _zero_ asymptotically as \(\Delta t\) increases and that this trend is anyway similar to that of \(exp(-\lambda\Delta t)\), even though it is not exactly equal to it. The last observation point with this method remains undefined, however, as it always yields \(log(0)\)A further clarification is that I avoided any consideration relative to data uncertainties, as this would distract us from examining the data trends. It is clear, however, that, as \(\Delta t\) increases, the uncertainties become very large. In Figure 1, both the normalized survival probability \(SP(\Delta t)\) (left) and the normalized \(\lambda(\Delta t)\) (right) are shown for DS1-DS4. A dashed line indicates the corresponding \(SP_{Poiss}(\Delta t)\). It is evident that all the survival functions do not resemble a Poissonian decay, even the de-clustered \(M\geq 8.50\) one. In this case, if we focus on \(\lambda(\Delta t)\), we observe that it is actually (quite notably) increasing, as the graph is in logarithmic scale along both axes. However, the de-clustered pattern is totally different from that of those that were not de-clustered. As stated before, Reference [11] showed that it is impossible to distinguish between the background seismicity, main events, and aftershocks. I go a step further, showing that the future event occurrence probability is strongly influenced by the history of the time series, and that clustered events actually lead to a totally different occurrence probability at a long temporal distance from the cluster. Figure 1: Survival probability (**left**) and \(\lambda(\Delta t)\) (**right**) for datasets 1–4 (DS1–DS4). Solid lines refer to empirical functions, while dashed lines refer to Poissonian ones. First of all, we need to consider that, if a dataset is non-Poissonian, it must have some memory in itself; this means that the classical representation, in terms of the survival probability function (as shown in Figure 1, left column), does not contain all the information required to describe the system's evolution. To address this issue, let us concentrate on DS1 and DS4, the only very large datasets available. Figures 2 and 3 show the distribution between earthquake triplets; in other words, for every \(\Delta t_{[i,i-1]}\) plotted along the abscissa, \(\Delta t_{[i+1,i]}\) is plotted along the ordinate.