Glitch time series and size distributions in eight prolific pulsars

TL;DR

This study analyzes glitch size and timing distributions in eight pulsars, revealing unique correlations in two pulsars and diverse statistical behaviors, advancing understanding of pulsar glitch mechanisms.

Contribution

It provides the first detailed statistical modeling of glitch properties across multiple pulsars, highlighting the unique correlation in PSR J0537-6910 and testing hypotheses about glitch populations.

Findings

Vela and PSR J0537-6910 show clear preferred glitch sizes and intervals.

PSR J0537-6910 uniquely exhibits a strong positive correlation between glitch size and waiting time.

Most pulsars have exponential waiting time distributions and diverse glitch size distributions.

Abstract

Glitches are rare spin-up events that punctuate the smooth slow-down of the rotation of pulsars. For the Vela pulsar and PSR J0537-6910, the glitch sizes and the times between consecutive events have clear preferred scales, contrary to the handful of other pulsars with enough glitches for such a study. Moreover, PSR J0537-6910 is the only pulsar showing a strong positive correlation between the size of each glitch and the waiting time until the following one. We attempt to understand this behavior through a detailed study of the distributions and correlations of glitch properties for the eight pulsars with at least ten detected glitches. We model the distributions of glitch sizes and times between consecutive glitches for this sample and perform Monte Carlo simulations to explore two hypotheses that could explain why the correlation is so much weaker in other pulsars than in PSR J0537-6910. We confirm the above results for the Vela pulsar and PSR J0537-6910, and verify that the latter is the only pulsar with a strong correlation between glitch size and waiting time to the following glitch. For the remaining six pulsars, the waiting time distributions are best fitted by exponentials, and the size distributions either by power laws, exponentials, or log-normal functions. Some pulsars in the sample yield significant Pearson and Spearman coefficients ($r_p$ and $r_s$) for the aforementioned correlation. Moreover, for all except the Crab, both coefficients are positive. Our simulations show that the weaker correlations in pulsars other than PSR J0537$-$6910 cannot be due to missing glitches too small to be detected. We also tested the hypothesis that each pulsar may have two kinds of glitches, namely large, correlated ones and small, uncorrelated ones. We found for all pulsars that their glitches are not consistent with a pure uncorrelated population either.

Tables

Table 1. Table 1 : Distributions of glitch sizes: results of the fits and the AIC weights for each model; using glitches with Δ ​ ν ≥ 0.01 ​ μ ​ Hz Δ 𝜈 0.01 𝜇 Hz \Delta\nu\geq 0.01\,\rm{\mu Hz} .

PSR Name	$w^{\mathrm{Gauss}}$	$w^{\text{Power law}}$	$w^{\text{L-N}}$	$w^{\mathrm{Exp}}$	$\hat{\mu }$	$\hat{\sigma }$	$\hat{\alpha }$	${\hat{\mu }}_{L-N}$	${\hat{\sigma }}_{L-N}$	$\hat{\lambda }$
					$\mu \mathrm{Hz}$	$\mu \mathrm{Hz}$				${\left(\mu \mathrm{Hz}\right)}^{-1}$
J0205$+$6449	${10}^{-8}$	$0.66$	0.33	${10}^{-5}$	15(5)	20(4)	1.27(6)	0.7(7)	2.5(3)	0.07(6)
B0531$+$21	${10}^{-17}$	0.02	$0.97$	${10}^{-7}$	1.2(5)	3(1)	1.4(1)	-1.3(3)	1.5(2)	0.8(7)
J0537$-$6910	$0.96$	${10}^{-24}$	${10}^{-8}$	0.03	15(1)	9.9(9)	1.19(5)	2.2(2)	1.3(2)	0.063(6)
J0631$+$1036	${10}^{-12}$	$0.94$	0.05	${10}^{-8}$	1(1)	3(1)	1.4(1)	-1.9(6)	2.1(4)	0.61(4)
B0833$-$45	$0.997$	${10}^{-13}$	${10}^{-6}$	0.002	21(2)	9(1)	1.2(4)	2.7(2)	1.2(4)	0.05(1)
B1338$-$62	${10}^{-5}$	0.07	$0.53$	0.4	2.5(5)	2.7(3)	1.36(5)	-0.1(3)	1.6(1)	0.4(1)
B1737$-$30	${10}^{-14}$	$0.82$	0.17	${10}^{-7}$	0.6(2)	1.0(2)	1.38(6)	-2.0(3)	1.9(1)	1.5(8)
B1758$-$23	0.06	0.004	0.07	$0.866$	0.6(1)	0.51(8)	1.3(2)	-1.2(4)	1.5(3)	1.7(6)

Table 2. Table 2 : Distributions of waiting times between successive glitches: results of the fits and the AIC weights for each model.

PSR Name	$w^{\mathrm{Gauss}}$	$w^{\mathrm{Powerlaw}}$	$w^{L-N}$	$w^{\mathrm{Exp}}$	$\hat{\mu }$	$\hat{\sigma }$	$\hat{\alpha }$	${\hat{\mu }}_{L-N}$	${\hat{\sigma }}_{L-N}$	$\hat{\lambda }$
					yr	yr				yr-1
J0205$+$6449	$0.001$	$0.40$	0.16	$0.43$	1.3(4)	1.4(4)	1.7(1)	-$0.2\left(3\right)$	1.0(1)	0.9(5)
B0531$+$21	${10}^{-4}$	${10}^{-5}$	0.15	$0.84$	1.3(2)	1.3(2)	1.4(1)	-$0.2\left(2\right)$	1.0(1)	0.8(2)
J0537$-$6910	$0.72$	${10}^{-10}$	0.07	0.2	0.28(2)	0.15(1)	1.64(8)	-1.44(9)	0.65(6)	4.3(4)
J0631$+$1036	${10}^{-4}$	${10}^{-5}$	0.20	$0.79$	1.4(4)	1.7(6)	1.3(2)	-0.3(3)	1.2(2)	0.7(3)
B0833$-$45	$0.993$	${10}^{-10}$	${10}^{-4}$	0.006	2.5(2)	1.2(1)	1.3(3)	0.7(2)	0.9(2)	0.41(9)
B1338$-$62	0.25	${10}^{-3}$	0.20	$0.54$	0.88(9)	0.42(4)	1.9(2)	-0.3(1)	0.51(5)	1.7(3)
B1737$-$30	${10}^{-5}$	${10}^{-6}$	0.17	$0.82$	0.9(1)	0.9(1)	1.44(7)	-0.6(1)	1.0(1)	1.2(2)
B1758$-$23	0.04	0.16	0.08	$0.72$	2.4(4)	1.4(2)	2.1(2)	0.7(1)	0.61(8)	0.6(2)

Table 3. Table 3 : Correlation coefficients between Δ ​ ν k Δ subscript 𝜈 𝑘 \Delta\nu_{k} and Δ ​ τ k + 1 Δ subscript τ 𝑘 1 \Delta\uptau_{k+1} .

PSR Name	$N_g$	$r_p$	$p_p$	$r_s$	$p_s$
J0205$+$6449	13	0.88	0.0002	$0.76$	0.004
B0531$+$21	25	-0.10	0.62	-0.12	0.57
J0537$-$6910	45	0.95	${10}^{-22}$	0.95	${10}^{-23}$
J0631$+$1036	17	0.93	${10}^{-7}$	0.20	0.45
B0833$-$45	20	0.24	0.31	0.31	0.21
B1338$-$62	23	0.59	0.003	0.70	0.0002
B1737$-$30	36	0.29	0.09	0.29	0.08
B1758$-$23	13	0.76	0.003	0.80	0.001

Table 4. Table 4 : Correlation coefficients for the pairs of variables ( Δ ​ ν k , Δ ​ τ k ) Δ subscript 𝜈 𝑘 Δ subscript τ 𝑘 (\Delta\nu_{k},\,\Delta\uptau_{k}) , and ( Δ ​ ν k , Δ ​ ν k − 1 ) Δ subscript 𝜈 𝑘 Δ subscript 𝜈 𝑘 1 (\Delta\nu_{k},\,\Delta\nu_{k-1}) .

	$\Delta {\nu }_k$ vs $\Delta {\tau }_k$
PSR Name	$r_p$	$p_p$	$r_s$	$p_s$
J0205$+$6449	$0.16$	$0.60$	$0.44$	$0.15$
B0531$+$21	$-0.02$	$0.90$	$0.40$	$0.05$
J0537$-$6910	$-0.08$	$0.60$	$-0.12$	$0.41$
J0631$+$1036	$-0.10$	$0.68$	$-0.18$	$0.49$
B0833$-$45	$0.55$	$0.01$	$0.27$	$0.24$
B1338$-$62	$-0.30$	$0.16$	$-0.18$	$0.41$
B1737$-$30	$-0.02$	$0.89$	$-0.10$	$0.56$
B1758$-$23	$-0.02$	$0.94$	$-0.04$	$0.89$

	$\Delta {\nu }_k$ vs $\Delta {\nu }_{k-1}$
PSR Name	$r_p$	$p_p$	$r_s$	$p_s$
J0205$+$6449	$-0.06$	$0.83$	$0.25$	$0.42$
B0531$+$21	$-0.10$	$0.61$	$-0.15$	$0.47$
J0537$-$6910	$-0.13$	$0.38$	$-0.16$	$0.29$
J0631$+$1036	$-0.12$	$0.65$	$0.32$	$0.21$
B0833$-$45	$-0.08$	$0.71$	$-0.12$	$0.59$
B1338$-$62	$-0.33$	$0.13$	$-0.13$	$0.55$
B1737$-$30	$-0.11$	$0.50$	$0.03$	$0.85$
B1758$-$23	$-0.02$	$0.92$	$-0.04$	$0.89$

Equations

$$M(x|\mu ,\sigma )=C_{\mathrm{Gauss}}\exp \left[\frac{-(x-\mu {)}^2}{2{\sigma }^2}\right]\text{,}$$

$$M(x|\alpha )=\frac{\alpha -1}{x_{\min }}{\left(\frac{x}{x_{\min }}\right)}^{-\alpha }\text{,}$$

$$M(x|{\mu }_{\mathrm{L}-\mathrm{N}},{\sigma }_{\mathrm{L}-\mathrm{N}})=\frac{C_{\mathrm{L}-\mathrm{N}}}{x}\exp \left[\frac{-(\ln x-{\mu }_{\mathrm{L}-\mathrm{N}}{)}^2}{2{\sigma }_{\mathrm{L}-\mathrm{N}}^2}\right]\text{,}$$

$$M(x|\lambda )=\lambda \exp \left[-\lambda (x-x_{\min })\right]\text{.}$$

$$P(\ge 6|7)=P(6|7)+P(7|7)=\frac{1}{16}=0.0625.$$

$$\Delta {\text{uptau}}_{k+1}=C\Delta {\nu }_k.$$

$$\Delta {\text{uptau}}_{k+1}=\Delta {\nu }_k\cdot 10^x,$$

Full text

*Aims. * We attempt to understand this behaviour through a detailed study of the distributions and correlations of glitch properties for the eight pulsars with at least ten detected glitches.

*Methods. * We model the distributions of glitch sizes and times between consecutive glitches for this sample. Monte Carlo simulations are used to explore two hypotheses that could explain why the correlation is so much weaker in other pulsars than in PSR J0537$-$6910.

*Results. * We confirm the above results for the Vela pulsar and PSR J0537$-$6910, and verify that the latter is the only pulsar with a strong correlation between glitch size and waiting time to the following glitch. For the remaining six pulsars, the waiting time distributions are best fitted by exponentials, and the size distributions either by power laws, exponentials, or log-normal functions. Some pulsars in the sample yield significant Pearson and Spearman coefficients ($r_{p}$ and $r_{s}$) for the aforementioned correlation. Moreover, for all except the Crab, both coefficients are positive. For each coefficient taken separately, the probability of this happening by chance is $1/16$. Our simulations show that the weaker correlations in pulsars other than PSR J0537$-$6910 cannot be due to missing glitches too small to be detected. We also tested the hypothesis that each pulsar may have two kinds of glitches, namely large, correlated ones and small, uncorrelated ones. The best results are obtained for the Vela pulsar, which exhibits a correlation with $r_{p}=0.68$ ($p$-value$=0.003$) if its 2 smallest glitches are removed. The other pulsars are harder to accommodate under this hypothesis, but their glitches are not consistent with a pure uncorrelated population either. We also found that all pulsars in our sample, except the Crab, are consistent with the previously found constant ratio between glitch activity and spin-down rate, $\dot{\nu}_{g}/|\dot{\nu}|=0.010\pm 0.001$, even though some of them have not shown any large glitches.

*Conclusions. *

To explain these results, we speculate that, except in the case of the Crab pulsar, all glitches draw their angular momentum from a common reservoir (presumably a neutron superfluid component containing $\approx 1\%$ of the star’s moment of inertia), but two different trigger mechanisms could be active, a more deterministic one for larger glitches and a more random one for smaller ones.

1 Introduction

The rotation frequencies $\nu$ of pulsars generally decrease slowly in time, but occasionally experience sudden increases $\Delta\nu$ that are usually accompanied by increases in the absolute value of their spin-down rates, $\dot{\nu}$ (Radhakrishnan & Manchester 1969; Reichley & Downs 1969; Shemar & Lyne 1996). These spin-up events, known as glitches, are infrequent, not periodic, and cover a wide range of sizes (from $\Delta\nu/\nu\sim 10^{-11}$ to $\Delta\nu/\nu\sim 10^{-5}$; Espinoza et al. 2011; Yu et al. 2013). The mechanism that generates these events is not completely understood, but they are believed to be caused by angular momentum transfer from an internal neutron superfluid to the rest of the neutron star (Anderson & Itoh 1975).

Thanks to the few long-term monitoring campaigns that keep operating, some since the 1970s (e.g. Hobbs et al. 2004; Yu et al. 2013), the number of detected glitches has slowly increased, thereby improving the significance of statistical studies in pulsar populations. McKenna & Lyne (1990), Lyne et al. (2000), and Espinoza et al. (2011) showed that the glitch activity $\dot{\nu}_{\rm{g}}$ (defined as the mean frequency increment per unit of time due to glitches) correlates linearly with $|\dot{\nu}|$. They also found that young pulsars (using the characteristic age, $\tau_{c}=-\nu/2\dot{\nu}$, as a proxy for age), which also have the highest $|\dot{\nu}|$, exhibit glitches more often than older pulsars, with rates varying from about one glitch per year to one per decade among the young pulsars. Using a larger and unbiased sample, Fuentes et al. (2017) confirmed that the size distribution of all glitches in a large and representative sample of pulsars is multi-modal (recently also seen by Konar & Arjunwadkar 2014; Ashton et al. 2017), with at least two well-defined classes of glitches: large glitches in a relatively narrow range $\Delta\nu\sim(10-30)\,\rm{\mu Hz}$, and small glitches with a much wider distribution, from $\sim 10\,\mathrm{\mu Hz}$ down to at least $10^{-4}\,\mathrm{\mu Hz}$. Further, Fuentes et al. (2017) found that a constant ratio $\dot{\nu}_{\rm{g}}/|\dot{\nu}|=0.010\pm 0.001$ is consistent with the behaviour of nearly all rotation-powered pulsars and magnetars. The only exception are the (few) very young pulsars, which have the highest spin-down rates, such as the Crab pulsar (PSR B0531$+$21) and PSR B0540$-$69.

Because glitches are rare events, the number of known glitches in the vast majority of pulsars is not enough to perform robust statistical analyses on individual bases. This has made people focus on the few objects that have the largest numbers of detected glitches (about 10 pulsars). The statistical distributions of glitch sizes and times between consecutive glitches (waiting times), for the nine pulsars with more than five known glitches at the time, were studied by Melatos et al. (2008). They found that seven out of the nine pulsars exhibited power-law-like size distributions and exponential waiting time distributions. The distributions of the other two (PSRs J0537$-$6910 and B0833$-$45, the Vela pulsar) were better described by Gaussian functions, setting preferred sizes and time scales. These results have been further confirmed by Fulgenzi et al. (2017) and Howitt et al. (2018), who also found that there are at least two main behaviours among the glitching pulsars.

Correlations between glitch sizes and the times to the nearest glitches, either backward or forward, are naturally expected. We know that glitch activity is driven by the spin-down rate (Fuentes et al. 2017), which suggests that glitches are the release of some stress that builds up at a rate determined by $|\dot{\nu}|$. If the stress is completely released at each glitch, then one should expect a correlation between size and the time since the last glitch. Conversely, if glitches occur when a certain critical state is reached, one should expect a correlation between size and the time to the next glitch, as longer times would be needed to come back to the critical state after the largest glitches. Moreover, if both assumptions were indeed correct, glitches would all be of equal sizes and occur periodically. However, with the exception of PSR J0537$-$6910 (see below), no other pulsars have shown significant correlations between glitch sizes and the times to the nearest events (e.g. Wang et al. 2000; Yuan et al. 2010; Melatos et al. 2018). This may be partly due to small-number statistics and might improve in the future, provided a substantial number of pulsars continue to be monitored for glitches.

The case of PSR J0537$-$6910, however, is very clear. With more than 40 glitches detected in $\sim 13$ yr, the statistical conclusions about its behaviour are much more significant than for any other pulsar. As first reported by Middleditch et al. (2006), its glitch sizes exhibit a strong correlation with the waiting time to the following glitch (see also Antonopoulou et al. 2018; Ferdman et al. 2018, who confirmed the correlation using twice as much data).

Antonopoulou et al. (2018) interpret this behaviour as an indication that glitches in this pulsar occur only once some threshold is reached. Moreover, this behaviour would imply that not necessarily all the stress is released in the glitches, thereby giving rise to the variety of (unpredictable) glitch sizes observed and the lack of backward time correlation.

In this work we study the sequence of glitches in the pulsars with at least ten detected events, by characterizing their distributions of glitch sizes and waiting times between successive glitches. Also, we test two hypotheses to explain why most pulsars do not show a correlation between glitch size and time to the following glitch: the effects of undetected small glitches and the possibility that two different classes of glitches are present in each pulsar.

2 Pulsars with at least ten detected glitches

To date, there are eight pulsars with at least 10 detected glitches (Fig. 2). PSRs J0205$+$6449, B0531$+$21 (the Crab pulsar), B1737$-$30, B1758$-$23, and J0631$+$1036 have been observed regularly by the Jodrell Bank Observatory (JBO, Hobbs et al. 2004). PSR B1338$-$62 has been observed by the Parkes telescope, and the Vela pulsar has been observed by several telescopes, including Parkes, the Jet Propulsion Laboratory, and others in Australia and Southafrica (e.g. Downs 1981; McCulloch et al. 1987; Yu et al. 2013; Buchner 2013). PSR J0537$-$6910 is the only object in our sample not detected in the radio band and was observed for 13 years by the Rossi X-ray Timing Explorer (RXTE, Antonopoulou et al. 2018; Ferdman et al. 2018). Glitch epochs and sizes were taken from the JBO online glitch catalog 333http://www.jb.man.ac.uk/pulsar/glitches/gTable.html, where more information and the appropriate references for each measurement can be found.

Figures 2 and 3 show that the Vela pulsar and PSR J0537$-$6910 produce glitches of similar sizes, particularly large glitches ($\Delta\nu>10$ $\mu$Hz), and in fairly regular time intervals. The absence of smaller glitches in these pulsars is not a selection effect, as it is quite unlikely that a considerable amount of glitches with sizes up to $\Delta\nu\sim 10\,\rm{\mu Hz}$, far above the detection limits reported in the literature (see Watts et al. 2015, and text below), could have gone undetected. On the other hand, the rest of the pulsars exhibit irregular waiting times and cover a wider range of sizes ($\Delta\nu\sim 10^{-3}-10$ $\mu$Hz).

The cadence of the timing observations varies considerably from pulsar to pulsar (and even with time for individual pulsars), and the sensitivity of the observations, from which the glitch measurements were performed, are also different between different pulsars. This means that the chances of detecting very small glitches are different for each pulsar and that the completeness of the samples towards small events might also be different (Espinoza et al. 2014). Nonetheless, in this study we use a single value to represent the glitch size below which samples are likely to be incomplete due to detectability issues. For an observing cadence of 30 days and a rotational noise of 0.01 rotational phases, glitch detection is severely compromised below sizes $\Delta\nu\sim 10^{-2}\,\rm{\mu Hz}$, especially if their frequency derivative steps are larger than $|\Delta\dot{\nu}|\sim 10^{-15}\,\rm{Hz\,s^{-1}}$ (see Watts et al. 2015). We use the above numbers to characterize the glitch detection capabilities in this sample of pulsars, but we note that such cadence and rotational noise are rather pessimistic values in some cases.

3 Distributions of glitch sizes and times between glitches

In the following, we model the distributions of glitch sizes ($\Delta\nu$, measured in $\mu$Hz) and the distributions of times between successive glitches ($\Delta\tau$, measured in yr) for each pulsar in our sample. Four probability density distributions are considered: Gaussian,

[TABLE]

power-law,

[TABLE]

log-normal,

[TABLE]

and exponential,

[TABLE]

The set $\{\mu,\sigma,\alpha,\mu_{\rm{L-N}},\sigma_{\rm{L-N}},\lambda\}$ are the fitting parameters. All the distributions are normalized in the range $x_{\rm{min}}$ to $\infty$. Formally, $x_{\rm{min}}$ is given by detection limits. However, it is not simple to define precise values for $\Delta\nu_{\rm{min}}$ and $\Delta\uptau_{\rm{min}}$ for each pulsar. Thus we use $\Delta\nu_{\rm{min}}=10^{-2}\,\mu$Hz for the glitch sizes (see previous section), and the smallest interval of time between glitches in each pulsar as $\Delta\uptau_{\rm{min}}$.

For the Gaussian and log-normal distributions the normalization constants $C_{\rm{Gauss}}$ and $C_{\rm{L-N}}$ were found numerically. We use the maximum likelihood technique to obtain the parameters of the models that describe best the data, and use the Akaike Information Criterion (AIC, Akaike 1974) to compare the different models (see also the Appendix in Fuentes et al. 2017).

Figures 4-5 and Tables 1-2 summarize the results of fitting these distributions to each pulsar. There is no single distribution type that can simultaneously describe all the pulsars satisfactorily, for either sizes or waiting times. The size distributions present a large variety (as also found in the model of Carlin & Melatos 2019): the log-normal distribution gives the best fit for the Crab pulsar and PSR B1338$-$62, power-law for PSRs J0631$+$1036, B1737$-$30, and J0205$+$6449, and exponential for PSRs B1758$-$23.

We also note that PSR J0205+6449 and PSR B1758$-$23 are the pulsars with the fewest recorded glitches in the sample (both have 13 glitches detected), hence we ought to wait and confirm this result once more events are detected.

In the case of PSRs J0537$-$6910 and B0833$-$45 (Vela), the best fit for both size and waiting time distributions are Gaussian functions. Their size distributions are centered at large sizes $\Delta\nu\approx 15$ and $20\,\rm{\mu Hz}$, respectively, consistent with the peak of large glitches in the combined distribution for all pulsars (Fuentes et al. 2017).

The distributions of times between successive glitches offer more homogeneous results. Besides the case of PSR J0537$-$6910 and the Vela pulsar (best modelled by Gaussian functions), the waiting time distributions for all the other pulsars are best represented by exponential functions. These results are in agreement with Melatos et al. (2008); Wang et al. (2012), and Howitt et al. (2018) for almost all the pulsars studied. The only exception is PSR B1338$-$62, for which Howitt et al. (2018) reported a local maximum in the distribution and classified this pulsar as a quasi-periodic glitcher.

If $\Delta\nu_{\rm{min}}$ is set to the size of the smallest detected glitch in each pulsar (rather than to $10^{-2}\,\mu$Hz), the results of the fits are very similar, and give parameters within the uncertainties presented in Table 1.

4 Time series correlations: Glitch size and time to the next glitch

Different studies have shown that for PSR J0537$-$6910 the glitch magnitudes $\Delta\nu_{k}$ are strongly correlated with the waiting times to the following glitch $\Delta\uptau_{k+1}$ (Middleditch et al. 2006; Antonopoulou et al. 2018; Ferdman et al. 2018, and see Fig. 6).

We test whether this correlation is also present in the other pulsars of the sample, and show the results in Table 6 and Fig. 6 (this is fairly consistent with Melatos et al. 2018, though we note that the samples of glitches are not exactly the same). None of them exhibits a correlation as clear as PSR J0537$-$6910. However, for PSRs J0205$+$6449, J0631$+$1036, B1338$-$62, and B1758$-$23, the Pearson correlation coefficients are larger than $0.5$ and the $p$-values are $\sim 10^{-3}$, or less. Therefore, at $95\%$ confidence level ($p$-values $<0.05$), we can reject the null hypothesis that $\Delta\nu_{k}$ and $\Delta\uptau_{k+1}$ are uncorrelated in these pulsars. Since the Pearson coefficient can be dominated by outliers, we also compute the Spearman rank correlation coefficient, obtaining similar or even stronger correlations, except for PSR J0631$+$1036.

It is also interesting to note that not only for PSR J0537$-$6910, but for all pulsars in the sample except the Crab, both the Pearson and Spearman correlation coefficients are positive. The probability of finding at least six out of seven pulsars having the same sign as our reference case, just by chance, is rather low. The probability of getting exactly $k$ successes among $n$ trials, with $1/2$ success probability in each trial, is $P(k\,|\,n)={n\choose k}(1/2)^{n}$. Thus, the probability of getting at least 6 successes in 7 trials is

[TABLE]

This low probability suggests that the waiting time to the following glitch is at least partially regulated by the size of the previous glitch.

In order to explain why the correlation for all other pulsars is much less clear than for PSR J0537$-$6910, we explore two hypotheses, both of which are motivated by noting that most glitches in PSR J0537$-$6910 are large:

• (I)

The correlation is intrinsically present in the full population of glitches of each pulsar, but glitches below a certain size threshold are not detected, thereby increasing by random amounts the times between the detected ones and worsening the correlation.

• (II)

There are two classes of glitches: glitches above a certain threshold size that follow the correlation, and glitches below the same threshold that are uncorrelated.

4.1 Hypothesis I: Incompleteness of the sample

In order to test the first hypothesis, we simulate a hypothetical pulsar with 100 glitches that follow a perfect correlation between $\Delta\nu_{k}$ and $\Delta\uptau_{k+1}$. The events smaller than a certain value are then removed to understand the effect of their absence in the correlation. The procedure is the following:

Glitch sizes are generated from a power-law distribution given by $dN/d\Delta\nu\propto\Delta\nu^{-\alpha}$, with power-law index $\alpha>1$. We choose a power-law distribution because it mainly produces small events, and we want to see the effect of removing a substantial fraction of them. Several different choices for $\alpha$ were considered. Here we only show the results for $\alpha=1.2$ and 1.4, as they generate distributions that resemble some of the ones observed.

The distributions do not have an upper cutoff, and the lower limit was varied so that, after reducing the sample of glitches (as we explain in step 3 below), the resulting sample covers the typical observed range of glitch sizes ($10^{-2}-10^{2}\,\rm{\mu Hz}$). 2. 2.

The time to the next glitch $\Delta\uptau_{k+1}$ is computed in terms of the glitch size $\Delta\nu_{k}$ as:

[TABLE]

The value of the proportionality constant $C$ is irrelevant in this case, since we are simulating a generic pulsar. 3. 3.

Steps (1) and (2) are repeated until a sequence of 100 glitches is reached. Then the 80 smallest are removed, thereby leaving a reduced sample of 20 to be analyzed, which is comparable to the number of glitches observed in each of our 8 pulsars. The lower limit for the distribution is computed analytically so that, after reducing the sample of glitches, the final sample covers the typical observed range of glitch sizes ($10^{-2}-10^{2}\,\rm{\mu Hz}$). 4. 4.

Finally, we calculate the time interval between each pair of successive glitches in the reduced sample, and determine both the Spearman and Pearson correlation coefficients between $\Delta\nu_{k}$ and $\Delta\uptau_{k+1}$.

After simulating $10^{4}$ cases, it was found that removing all glitches smaller than a certain value has a minor effect on the correlation. Representative realizations are shown in Fig. 7, where the correlation between $\Delta\nu_{k}$ and $\Delta\uptau_{k+1}$ is plotted in log-scale to show more clearly the dispersion produced by the removal of the smallest glitches. We observe that missing small glitches does not substantially worsen the correlation: more than $90\%$ of the realizations give correlation coefficients $\geq 0.95$ (both Pearson and Spearman).

For $\alpha>1.4$ the distribution becomes narrower, accumulating towards the lower limit. Since a large fraction of the simulated glitches have very similar sizes, after removing the 80 smallest glitches the correlation does worsen, and yields correlation coefficients between $0.4$ and $0.9$, which are similar to those exhibited by the real data. However, in these cases the distributions of glitch sizes differ strongly from those observed for the pulsars in our sample.

From these simulations, we conclude that it is unlikely that the non-detection of all the glitches below a certain detection limit is the explanation for the low observed correlations in pulsars other than PSR J0537$-$6910.

4.2 Hypothesis II: Two classes of intrinsically different glitches

The second hypothesis states that pulsars exhibit two classes of glitches: larger events, which follow a linear correlation between $\Delta\nu_{k}$ and $\Delta\uptau_{k+1}$; and smaller events, for which these variables are uncorrelated. We allow the point of separation between large and small glitches to be different for each pulsar.

To visualize whether this hypothesis works, correlation coefficients (for the same pair of variables, $\Delta\nu_{k}$ and $\Delta\uptau_{k+1}$) were calculated for sub-sets of glitches of the original sample. The sub-sets are defined as all glitches with sizes larger or equal to a given $\Delta\nu_{\textrm{min}}$. Correlation coefficients as a function of $\Delta\nu_{\textrm{min}}$ are plotted in Fig. 8 for each pulsar. Visual inspection of the plots immediately tells us that by removing small glitches no pulsar reaches the level of correlation observed for PSR J0537$-$6910, for both correlation tests.

In the following we explore the curves in Fig. 8 in some more detail. For that purpose, Monte Carlo simulations of pulsars with correlated and uncorrelated glitches were performed. Since the underlying glitch size distributions of the pulsars in the sample are unknown, we use the measured values of a given pulsar. The following is the procedure for one realization:

The glitches larger than a certain value $\Delta\nu^{\star}$ are chosen in random order and assigned epochs according to their size. The first one is set at an arbitrary epoch and the epochs of the following ones are assigned according to

[TABLE]

where $x$ is drawn from a Gaussian distribution centred at $\bar{x}=\log(C)$ and with a standard deviation equal to $\sigma_{\bar{x}}$. The latter allows us to introduce a dispersion in the correlation of the simulated glitches. The distribution of $\log(\Delta\uptau_{k+1}/\Delta\nu_{k})$ for all glitches with $\Delta\nu>5\,\mu$Hz in PSR J0537$-$6910 can be well modelled by a Gaussian distribution with standard deviation $\sigma_{0537}=0.085$ (in logarithmic scale, if $\Delta\uptau_{k+1}$ is measured in days and $\Delta\nu_{k}$ is measured in $\mu$Hz). In the simulations, $\sigma_{\bar{x}}$ was set either to zero (i.e. $x=\log(C)$, perfect correlation) or to multiples of $\sigma_{0537}$. 2. 2.

The glitches smaller than $\Delta\nu^{\star}$ are distributed randomly over the time span between the first and the last correlated glitches. The resulting waiting times of all, correlated and uncorrelated glitches are then multiplied by a factor that ensures that their sum equals the time in between the first and the last observed glitches. 3. 3.

Steps 1 and 2 were repeated $10^{4}$ times for each considered value of $\Delta\nu^{\star}$.

The plots in Fig. 9 show the results of simulations using the glitch sizes of PSR J0537$-$6910 and $\sigma_{\bar{x}}=\sigma_{0537}$ for three values of $\Delta\nu^{\star}$. The results are shown via curves of $r$ versus $\Delta\nu_{\textrm{min}}$, to compare with Fig. 8. The shaded areas represent the $70\%$ of the correlation coefficients closer to the median of all realizations. We visually inspected the distributions of $r_{p}$ and $r_{s}$ for all possible $\Delta\nu_{\textrm{min}}$ values, and for many $\Delta\nu^{\star}$ cases. It was verified that the median is sufficiently close to the maximum of the distribution in most cases. Though, this tends to fail for the largest $\Delta\nu_{\textrm{min}}$ values, where the $r_{p}$ and $r_{s}$ distributions are rather flat. But this is irrelevant because any conclusion pointing to a case in which only a few glitches are correlated (large $\Delta\nu_{\textrm{min}}$) would have little statistical value, regardless of the above. Thus we are confident that the shaded areas effectively cover the most possible outcomes of series of glitches under the assumptions considered.

We now use the plots in Fig. 9 to understand the curves of the correlation coefficients as functions of $\Delta\nu_{\textrm{min}}$ in Fig. 8, in the frame of hypothesis II:

• (a)

If all glitches were correlated, which is the case shown in the leftmost plot in Fig. 9, the correlation coefficients would decrease gradually as $\Delta\nu_{\textrm{min}}$ increases. This is because a progressive reduction of the sample, starting from the smallest events (i.e. increasing the remaining waiting times by small random amounts), will gradually kill the correlation. Note that the correlation coefficients of the simulated glitches start at values just below $1.0$ for the smallest $\Delta\nu_{\textrm{min}}$, just like the observations of PSR J0537$-$6910. This is because $\sigma_{\bar{x}}=\sigma_{0537}$ in those simulations. Only for $\sigma_{\bar{x}}=0$ the simulations would start at correlation coefficients equal to $1.0$.

• (b)

If only glitches above a certain size $\Delta\nu^{\star}$ were correlated, the correlation coefficients would improve as small glitches are eliminated, and the remaining sub-set approaches the one in which all glitches are correlated (as in the middle plot of Fig. 9). One would expect a maximum correlation for $\Delta\nu_{\textrm{min}}\sim\Delta\nu^{\star}$, and a gradual decrease as $\Delta\nu_{\textrm{min}}$ increases beyond $\Delta\nu^{\star}$.

• (c)

If there were no correlated glitches, we should expect a rather flat curve of low correlation coefficients oscillating around zero (rightmost plot in Fig. 9).

The behaviours just described correspond to the general trends exhibited by the shaded areas in Fig. 9, which evolve smoothly with $\Delta\nu_{\textrm{min}}$. However, particular realizations show abrupt variations, of both signs, just as the observations do in Fig. 8.

Clearly, PSR J0537$-$6910 is best represented by case (a). Indeed, both correlation coefficients for this pulsar are maximum (and very similar) when all glitches are included and they decrease gradually as the smallest glitches are removed (Fig. 8). Nonetheless, we note that $r_{p}$ stays above $0.9$ (and $p_{p}<3\times 10^{-12}$) for $\Delta\nu_{\textrm{min}}\leq 7\,\mu$Hz, hence it is possible that the smallest glitches are not correlated. Another indication for this possibility is that the six glitches below $5\,\mu$Hz fall to the right of the distribution of $\log(\Delta\uptau_{k+1}/\Delta\nu_{k})$ for all glitches, and the width of the distribution is reduced considerably (from more than 2 decades to a half decade) when they are removed. In other words, the straight line that best fits the ($\Delta\uptau_{k+1}$, $\Delta\nu_{k}$) points passes closer to the origin (a more physically motivated situation, Antonopoulou et al. 2018), and the data exhibit a smaller dispersion around this line, when the smallest glitches are not included.

The pulsars B1338$-$62, and B1758$-$23 may in principle also correspond to case (a). As mentioned at the beginning of section 4, they present mildly significant correlations when all their glitches are considered, and both their $r_{p}$ and $r_{s}$ curves in Fig. 8 decrease as $\Delta\nu_{\textrm{min}}$ increases. By performing simulations with $\Delta\nu^{\star}=0$, and for different values of $\sigma_{\bar{x}}$, we find that the correlation coefficients of PSR B1758$-$23 are within the range of $70\%$ of the possible outcomes if $\sigma_{\bar{x}}$ is set to 5-6 times $\sigma_{0537}$.

For PSR B1338$-$62 the situation is less clear because the amplitudes of the variations of both $r_{p}$ and $r_{s}$ for $\Delta\nu_{\textrm{min}}<1\,\mu$Hz are rather high. One possible interpretation is that all glitches are correlated and the variations are due to the correlation not being perfect (i.e. $\sigma_{\bar{x}}\neq 0$). We find that only for $\sigma_{\bar{x}}\geq 10\times\sigma_{0537}$ the simulations can reproduce such behaviour and the observed values. Another possibility is that $\Delta\nu^{\star}\sim 0.2\,\mu$Hz, which could explain the local maxima of $r_{p}$ and $r_{s}$ around that value. The maxima and subsequent values can indeed be reproduced with lower levels of noise, $\sigma_{\bar{x}}=5\times\sigma_{0537}$. But for smaller values of $\Delta\nu_{\textrm{min}}$ most realizations ($>70\%$) give correlation coefficients below $0.5$, thus they fail at reproducing the observed $0.6$-$0.7$ at $\Delta\nu_{\textrm{min}}=0$.

It is clear that Hypothesis II does not apply to this pulsar directly, and that the observations are not consistent with a set of uncorrelated glitches either. Based on the lack of glitches with sizes equal or less than $0.1\,\mu$Hz after MJD $\sim$ 50400 (Fig. 2), we speculate that the sample might be incomplete for glitches smaller than this size after this date777This would be a more extreme case than those considered for the Hypothesis I because $0.1\,\mu$Hz is a rather high limit..

The pulsars J0205$+$6449 and J0631$+$1036 also exhibit significant Pearson correlations when all their glitches are considered. However, their $r_{s}$ curves tend to increase with $\Delta\nu_{\textrm{min}}$ rather to decrease. As mentioned before, the Pearson test can be affected by outliers, hence the behaviour we see for $r_{p}$ is likely due to the very broad size and waiting times distributions and the low numbers of events towards the high ends of the distributions, which produce outlier points for both pulsars (Fig.6). It is therefore difficult to conclude anything for PSR J0631+1036. Moreover, the observed behaviour is very hard to reproduce by the simulations, even for high levels of noise (we tried up to $\sigma_{\bar{x}}=12\times\sigma_{0537}$). Perhaps its largest glitches ($\Delta\nu\geq 0.1\,\mu$Hz) are indeed correlated, but the statistics are too low to conclude anything.

For J0205, however, the Spearman coefficients $r_{s}$ are rather high ($>0.55$ for all $\Delta\nu_{\textrm{min}}$) and both coefficients become similar and even higher for $\Delta\nu_{\textrm{min}}>1\,\mu$Hz. It is possible that glitches above this size are correlated in this pulsar. We find that the observed $r_{p}$ and $r_{s}$, and their evolution with $\Delta\nu_{\textrm{min}}$, are within the $70\%$ of simulations with $\Delta\nu^{\star}=1.3\,\mu$Hz and for $\sigma_{\bar{x}}=2\times\sigma_{0537}$. We note, however, that in this case the correlation coefficients observed for $\Delta\nu_{\textrm{min}}\leq 0.1\,\mu$Hz are higher than the vast majority of the realizations. Perhaps the small glitches are also correlated and follow their own relation, though we did not simulate such scenario. We conclude that the Hypothesis II does not fully explain this pulsar, although the 8 glitches above $1\,\mu$Hz appear to be well correlated indeed.

The Vela pulsar is the only pulsar in the sample that seems well represented by case (b). The highest $r_{p}=0.68$ has a probability $p_{p}=0.003$ and is obtained for $\Delta\nu_{\textrm{min}}\sim 2\,\mu$Hz. Both $r_{p}$ and $r_{s}$ decline monotonically for larger $\Delta\nu_{\textrm{min}}$ values. This behaviour suggests that glitches of sizes above $\sim 2\,\mu$Hz might indeed be correlated, but the correlation is somewhat noisy. The observed correlation coefficients fall within the middle 70$\%$ of the realizations if $\sigma_{\bar{x}}=2\times\sigma_{0537}$ and for $\Delta\nu^{\star}=2$-$10\,\mu$Hz. The case $\Delta\nu^{\star}=9.35\,\mu$Hz is presented in Fig. 10. We prefer this case because simulations for $\Delta\nu^{\star}=2\,\mu$Hz tend to fail at reproducing the low correlation coefficients ($\leq 0.4$) observed for the smallest $\Delta\nu_{\textrm{min}}$.

Finally, the cases of PSRs B0531+21 (the Crab) and B1737$-$30 are rather inconclusive. The Crab pulsar is perhaps the pulsar for which case (c) applies the best. Both correlation coefficients are negative or positive, and in both cases stay at relatively low absolute values, which leads to the conclusion that there are no correlated glitches in the Crab pulsar. We note that the high $r_{p}$ and $r_{s}$ values observed for $\Delta\nu_{\textrm{min}}\sim 0.6\,\mu$Hz are obtained with the 5-6 largest events and that a linear fit to their $\Delta\nu_{k}-\Delta\uptau_{k+1}$ does not pass close to the origin.

The case of B1737$-$30 is more complex. The observations show two $\Delta\nu_{\textrm{min}}$ values, $0.0015$ and $0.03\,\mu$Hz, after which the correlation coefficients decrease with the removal of more small glitches (Fig. 8). This behaviour is hard to reproduce under Hypothesis II, unless the dispersion of the correlation is increased considerably, to $10\times\sigma_{0537}$ or more. We conclude that Hypothesis II does not apply to this pulsar directly and that there is some extra complexity, as the data are also inconsistent with a set of purely uncorrelated glitches.

Surprisingly, even though no pulsar complies perfectly with Hypothesis II, and the only way in some cases is to increase the dispersion of the correlation ($\sigma_{\bar{x}}\gg\sigma_{0537}$), there is no pulsar in the sample that is well represented by case (c) (only the Crab, to some extent).

Therefore, the sizes of at least some glitches must be positively correlated with the times to the next glitch in the available datasets. The question is why this correlation is much stronger in PSR J0537$-$6910 than in all other pulsars of our sample. Could this be an effect of its particularly high spin-down rate? Or the fact that most of its glitches are large? It could be that the correlations are indeed there, as stated in Hypothesis II, but for some reason exhibit high $\sigma_{\bar{x}}$ values. Maybe the fact that the glitches in PSR J0537$-$6910 occur so frequently ensures that the relationship stays pure. But it could also be that reality was more complex. For instance, it could be that both small and large glitches were correlated, but each of them followed a different law.

5 Other correlations

We looked for other correlations between the glitch sizes and the times between them. Specifically, we tried $\Delta\nu_{k}$ vs $\Delta\uptau_{k}$ (size of the glitch versus the time since the preceding glitch), and $\Delta\nu_{k}$ vs $\Delta\nu_{k-1}$ (size of the glitch versus the size of the previous glitch). No pulsar shows a significant correlation between these quantities (Table 4).

Nearly all pulsars in our sample show negative correlation coefficients (both, Pearson and Spearman) for $\Delta\nu_{k}$ vs $\Delta\uptau_{k}$. The only exceptions are the Vela pulsar and PSR J0205$+$6449 (see Table 4). Our results are in general agreement with Melatos et al. (2018), who also found a lack of correlation between $\Delta\nu_{k}$ vs $\Delta\uptau_{k}$ for individual pulsars.

For $\Delta\nu_{k}$ vs $\Delta\nu_{k-1}$, in most cases the correlation coefficients are close to zero and the $p$-values are larger than $0.2$, i.e., no individual pulsar shows a significant correlation. However, the results could still be meaningful for the sample as a whole because all the pulsars have negative correlation coefficients, except for the Spearman coefficients for PSRs J0631$+$1036 and B1737$-$30). The probability of getting all Pearson’s correlations coefficients of the same sign just by chance, regardless of whether the sign is positive or negative, is $2\times p_{\mathrm{binom}}(8|8)=0.007$. This could establish an interesting constraint on the glitch mechanism: Smaller glitches are somewhat more likely to be followed by larger ones, and vice-versa. However, this statement has to be confirmed with more data in the future.

6 Discussion

Fuentes et al. (2017) found that all pulsars (with the strong exception of the Crab pulsar and PSR B0540$-$69) are consistent with a constant ratio between the glitch activity, $\dot{\nu}_{\rm g}$, and the spin-down rate, $\dot{\nu}_{\rm{g}}/|\dot{\nu}|=0.010\pm 0.001$, i.e., $\approx 1\%$ of their spin-down is recovered by the glitches. This fraction has been interpreted as the fraction of the moment of inertia in a superfluid component that transfers its angular momentum to the rest of the star in the glitches (Link et al. 1999; Andersson et al. 2012). Fuentes et al. (2017) used the observed bimodal distribution of glitch sizes to distinguish between large and small glitches, with the boundary at $\Delta\nu=10\,\mathrm{\mu Hz}$, and argued that the constant ratio is determined by the large glitches, whose rate, $\dot{N}_{\ell}$ is also proportional to $|\dot{\nu}|$. In this scenario, the much lower (sometimes null) glitch activities measured in many low-$|\dot{\nu}|$ pulsars are due to their observation time spans not being long enough to include any large glitches (or any glitch at all). Interestingly, the pulsars in our sample (except the Crab) are quite consistent with the constant ratio (Fig. 11), even those, like PSRs B1338$-$62, B1737$-$30, and B1758$-23$, which do not have any large glitches contributing to their activities.

On the other hand, pulsars with higher spin-down rates also have a larger fraction of large glitches. At the highest spin-down rates ($|\dot{\nu}|\geq 10^{-11}$ Hz s*-1*), the production of large glitches becomes comparable and sometimes higher than the production of small glitches, again with the notorious exception of the Crab and PSR B0540$-$69. This trend is also followed by the pulsars in our sample: all large glitches (but one in PSR J0631+1036), are concentrated in PSRs J0205$+$6449, J0537$-$6910, and the Vela pulsar, which are (together with the Crab) the ones with largest $|\dot{\nu}|$ values (see Fig. 2 and 11).

Thus, it seems to be the case that both large and small glitches draw from the same angular momentum reservoir (for all but the very young, Crab-like pulsars), but have different trigger mechanisms, the large ones being produced once a critical state is reached, whereas small ones occur in a more random fashion. For reasons still to be understood, the glitch activity of relatively younger, high $|\dot{\nu}|$, Vela-like pulsars is dominated by large glitches, whereas for smaller $|\dot{\nu}|$ the large glitches become less frequent, both in absolute terms and relative to the small ones (Wang et al. 2000; Espinoza et al. 2011). In this context, it is interesting to note that recent long-term braking index measurements indicate that Vela-like pulsars move towards the region where PSRs J0631+1036, B1737$-$30, and B1758$-$23 are located on the $P$–$\dot{P}$ diagram (Espinoza et al. 2017).

7 Summary and Conclusions

We studied the individual glitching behaviour of the eight pulsars that today have at least ten detected glitches. Our main conclusions are the following:

We confirm the previous result by Melatos et al. (2008) and Howitt et al. (2018) that, for Vela and PSR J0537$-$6910, the distributions of both their glitch sizes and waiting times are best fitted by Gaussians, indicating well-defined scales for both variables. For all other pulsars studied, the waiting time distribution is best fitted by an exponential (as would be expected for mutually uncorrelated events), but they have a variety of best-fitting size distributions: a power law for PSR J0205+6449, J0631+1036, and B1737$-$30, a log-normal for the Crab and PSR B1338$-$62, and an exponential for PSR B1758$-$23. 2. 2.

All pulsars in our sample, except for the Crab, have positive Spearman and Pearson correlation coefficients for the relation between the size of each glitch, $\Delta\nu_{k}$, and the waiting time to the following glitch, $\Delta\tau_{k+1}$. For each coefficient, the probability for this happening by chance is $1/16=6.25\%$. Both coefficients also stay positive as the small glitches are removed (see Fig. 8). 3. 3.

PSR J0537$-$6910 shows by far the strongest correlation between glitch size and waiting time until the following glitch ($r_{p}=r_{s}=0.95$, $p$-values $\lesssim 10^{-22}$). Another three pulsars, PSRs J0205$+$6449, B1338$-$62, and B1758$-$23, have quite significant correlations ($p$-values $\leq 0.004$ for both coefficients). 4. 4.

Our first hypothesis to explain the much weaker correlations in all other pulsars compared to PSR J0537$-$6910, namely missing glitches that are too small to be detected, is very unlikely to be correct. Our Monte Carlo simulations show that, for reasonable glitch size distributions, it cannot produce an effect as large as observed. 5. 5.

Our alternative hypothesis, namely that there are two classes of glitches, large correlated ones and small uncorrelated ones, comes closer to reproducing the observed relations; notably for PSRs J0205$+$6449 and Vela. The resulting correlations for both pulsars present dispersions that are twice the one observed for PSR J0537$-$6910. For the other pulsars, the required dispersions to accommodate this hypothesis are much larger. 6. 6.

The correlation coefficients between the sizes of two successive glitches, $\Delta\nu_{k-1}$ and $\Delta\nu_{k}$, as well as between the size of a glitch, $\Delta\nu_{k}$ and the waiting time since the previous glitch, $\Delta\tau_{k}$, are generally not significant in individual pulsars, but they are negative for most cases, suggesting some (weaker) relation also among these variables. 7. 7.

Except for the Crab, all pulsars in our sample are consistent with the constant ratio between glitch activity and spin-down rate, $\dot{\nu}_{\mathrm{g}}/|\dot{\nu}|=0.010\pm 0.001$ (Fuentes et al. 2017). This includes cases dominated by large glitches, as well as others with only small glitches. 8. 8.

The previous results suggest that large and small glitches draw their angular momentum from a common reservoir, although they might be triggered by different mechanisms. Large glitches, which dominate at large $|\dot{\nu}|$ (except for the Crab and PSR B0540$-$69), might occur once a certain critical state is reached, while small glitches, dominating in older pulsars with lower $|\dot{\nu}|$, occur at essentially random times.

All the above is based on the behaviour of the pulsars with the most detected glitches. Even though we have shown before that the activity of all pulsars appears to be consistent with one single trend, these pulsars could still be outliers among the general population. Only many more years of monitoring will clarify the universality of these results.

Acknowledgements.

We thank Vanessa Graber and Simon Guichandut for valuable comments on the first draft of this article. We are also grateful to Wilfredo Palma for conversations that guided us at the beginning of this work. We also thank Ben Shaw for information regarding the detection of recent glitches and for keeping the glitch catalog up to date. This work was supported in Chile by CONICYT, through the projects ALMA31140029, Basal AFB-170002, and FONDECYT/Regular 1171421 and 1150411. J.R.F. acknowledges partial support by an NSERC Discovery Grant awarded to A. Cumming at McGill University. C.M.E. acknowledges support by the Universidad de Santiago de Chile (USACH).