Anti-glitches in the Ultraluminous Accreting Pulsar NGC 300 ULX-1 Observed with NICER

Ultraluminous X-ray sources (ULXs) are extragalactic, bright X-ray point sources with super-Eddington luminosities for isotropic emission due to accretion onto a neutron star or stellar-mass black hole (see Kaaret et al. 2017, for a recent review). Following the discovery of X-ray pulsations from a handful of ULXs, with periods ranging from 0.42 s to ∼31 s (Bachetti et al. 2014; Fürst et al. 2016; Israel et al. 2017a, 2017b), it has been suggested that the central engines of most ULXs could be highly magnetized accreting neutron stars, based on theoretical (King & Lasota 2016; Mushtukov et al. 2017) or observational arguments (Koliopanos et al. 2017; Walton et al. 2018b). Further evidence comes from the claimed discovery of cyclotron resonant scattering features (CRSFs) in two ULXs (Brightman et al. 2018; Walton et al. 2018a), albeit one of them with no detected pulsations.

The ULX in the nearby spiral galaxy NGC 300 (D = 1.87 ± 0.12 Mpc; Rizzi et al. 2006) was mistakenly classified as a supernova (SN 2010da) upon its discovery (Monard 2010). Follow up Chandra observations showed that the system is instead a high-mass X-ray binary undergoing an outburst (Binder et al. 2011). Using observations obtained in 2016 December with XMM-Newton and NuSTAR, Carpano et al. (2018b) discovered X-ray pulsations at a frequency of ν = 31.6 mHz (P = 31.6 s, first reported in Carpano et al. 2018a), thus revealing the compact object as a neutron star. They also found that the pulsar is spinning up very rapidly with $\dot{\nu }\,=5.57\times {10}^{-10}$ Hz s⁻¹.

This article reports on the regular monitoring of NGC 300 ULX-1 with the Neutron Star Interior Composition Explorer (NICER; Gendreau et al. 2016), and the discovery of three apparent spin-down glitches that occurred over a span of about 120 days. In the many hundreds of radio, X-ray, and γ-ray pulsars monitored for glitches, all of the hundreds of detected glitches were spin-up glitches (Δν > 0; see, e.g., Espinoza et al. 2011; Yu et al. 2013; Fuentes et al. 2017). Among magnetars, observations reveal approximately twenty glitches, consisting of both spin-up and spin-down ones (see, e.g., Dib & Kaspi 2014). The largest magnetar spin-down glitch (with Δν = −21 μHz) is inferred to possibly have occurred during a giant flare of SGR 1900 + 14 (Woods et al. 1999), while the second largest was over a hundred times smaller (Δν ∼ −0.09 μHz) and accompanied by radiative changes (Archibald et al. 2013). The glitches described above all occur in nonaccreting pulsars, whose long-term spin frequencies decrease with time ($\dot{\nu }\lt 0$) due to rotational energy loss from electromagnetic dipole radiation (Pacini 1968; Gunn & Ostriker 1969). For pulsars that accrete matter from a binary companion, the torque due to accretion can cause the neutron star to spin up or spin down and is usually much stronger than the dipole radiation torque. There are only a few reported glitches in accreting pulsars; two spin-up glitches with Δν ≈ 1 μHz, each occurring during an accretion outburst episode (Galloway et al. 2004; Serim et al. 2017), and possibly two spin-down glitches (Klochkov et al. 2009).

In response to a series of Astronomer's Telegrams reporting the discovery of pulsations and the period evolution of NGC 300 ULX-1 (Bachetti et al. 2018a, 2018b; Carpano et al. 2018a; Grebenev & Mereminskiy 2018a, 2018b; Kennea 2018; Vasilopoulos et al. 2018a), NICER observed NGC 300 ULX-1 on 2018 February 6 for 1.6 ks. Following that, the source entered NICER's solar exclusion zone and became unobservable until 2018 May 1, when NGC 300 ULX-1 was again sufficiently far from the Sun. NICER then started a monitoring campaign on the source, primarily designed to discover evidence of an orbital period. NICER generally has a visibility window of a few hundred to ∼2000 s in one segment of an orbit, and observations were generally made during 1–10 consecutive orbital segments per day (most often 2, 3, or 4). Here we present an analysis of NICER ObsIDs 1034200101 (2018 February 6) through 1034200198 (2018 October 20).

The X-ray timing instrument of NICER provides ≈1900 cm² of collecting area at 1.5 keV thanks to its 56 coaligned X-ray concentrating optics paired with single-pixel silicon drift detectors (SDDs) working in the 0.2–12 keV energy band (52 detectors are operational on orbit). The GPS time tagging of photons permits obtaining an accuracy better than 100 ns root mean square (rms; LaMarr et al. 2016; Prigozhin et al. 2016).

The data were processed following standard procedures using the NICER data analysis software NICERDAS version 4, together with HEASOFT version 6.24. The standard filtering criteria define the good time intervals for which the NICER pointing offset is <0015 from the nominal source position; the source is at least 30 ° from the Earth limb (at least 40 ° in the case of a Sun-illuminated Earth), and intervals when NICER was not within the boundaries of the South Atlantic Anomaly. In a few observations, we masked "hot" detectors exhibiting out-of-family count rates. Also, for ObsID 1034200101, we had to remove MPU3¹⁹ because of extreme count rates caused by solar loading on the detectors, and relaxed the filtering for the bright Earth to 30 °. For each observation, we selected events in the 0.4–12 keV band and barycentered the data using the NASA/JPL DE405 solar system ephemeris and the coordinates of the ULX from Chandra imaging data (R.A. = 13770208, decl. = −37695417; Binder et al. 2011). All barycentric times and frequencies in this work use the Barycentric Dynamical Time (TDB) timescale. Some of the observations exhibited background flares when the International Space Station (ISS) was at extreme northern or southern latitudes. We applied a count rate cut of 7.5 c s⁻¹ on non-overlapping 8 s intervals to exclude the times of such flares.

The timing analysis demonstrates the power of NICER for timing studies of faint X-ray pulsars. NICER's flexible scheduling enabled monitoring of the source in a way that maximized the ability to do precise timing. We scheduled observations to cover two or more consecutive ISS orbits, spreading a relatively small amount of observing time over a duration long enough to make a precise frequency measurement, but not so long that the window function caused significant ambiguity in the pulse counting between data segments. These precise frequency measurements allowed us to bootstrap a fully coherent timing solution. In the following, we present the timing analysis of the pulsations, working from coarse frequency measurements to phase-coherent timing.

For the initial timing analysis, we binned the data at 1/128 s and used the PRESTO software package (Ransom et al. 2002) to perform a periodicity search over a range of possible accelerations. In essentially every observation until 2018 August 22, when the source seems to have faded significantly, the ULX pulsations (at around 53–57 mHz) were detected as the highest significance peak coming out of the search. This allowed us to confirm that the period evolution was continuing as previously reported. In addition, we confirmed the broad Gaussian pulse profile of this ULX reported previously (Carpano et al. 2018b). Figure 1 shows the pulse profile from MJD 58246, NICER's densest day of observations. The profile is reasonably well modeled as a Gaussian with a full width at half maximum (FWHM) of 0.35 cycles.

Zoom In Zoom Out Reset image size

To go from initial periodicity detections to precise period measurements, we employed a maximum likelihood technique to measure pulsed frequencies and significance. In this case, given a selection of photons from a segment of the observation, we searched over a grid of trial ν values. For each trial we computed the pulse phase of each photon according to a model with the trial ν (and $\dot{\nu }$ fixed at the long-term average, 4.16 × 10⁻¹⁰ Hz s⁻¹). We modeled the probability distribution of this set of phases using a Gaussian template with a FWHM 0.35 of the pulse period and subsequently maximized the likelihood to determine the optimal values of the Gaussian peak position and amplitude. By profiling over these nuisance parameters, we arrived at a profile likelihood for ν; see Kerr et al. (2015) and Ray et al. (2011) for a discussion of similar likelihood techniques. This finally yielded a maximum likelihood estimate for the pulse frequency at that observation, as well as an estimate of the measurement uncertainty determined by the frequency range over which the log likelihood decreased by 0.5.

The rapid and noisy period evolution, combined with the window function of the observations, made precision timing difficult, so we developed a bootstrapping technique to make precise pulse frequency measurements. The first step was making likelihood frequency measurements for each nearly contiguous segment of NICER data (i.e., with no gaps greater than 1000 s and a total span less than 3000 s). These measurements are low-precision, but unambiguous. We note that the last significant pulsed detection was at MJD 58352.5 (2018 August 22; ObsID 1034200183). All subsequent observations show no strong evidence of pulsations from NGC 300 ULX-1; this includes an 800 s observation at MJD 58351.5, followed by a gap of a month, then 21 observations between MJD 58386 and 58411.²⁰

The measured frequencies for each segment of NICER data are shown in Figure 2. Fitting the pulse frequency evolution to a simple model with two frequency derivatives gives the best fit $\nu ={\nu }_{0}+\dot{\nu }{(t-{t}_{0})+\tfrac{1}{2}\ddot{\nu }(t-{t}_{0})}^{2}$, where ν₀ = 53.350(5) mHz, $\dot{\nu }=4.16(5)\times {10}^{-10}$ Hz s⁻¹, $\ddot{\nu }=-2.5(12)\times {10}^{-18}$ Hz s⁻², and t₀ = MJD(TDB) 58243.515. Here, and elsewhere in the paper, the numbers in parentheses are the 1σ uncertainty in the last digit. If the spin down proceeds as $\dot{{\rm{\Omega }}}\propto {{\rm{\Omega }}}^{n}$ (where Ω is the angular velocity), the "braking index" n is defined as $n=\nu \ddot{\nu }/{\dot{\nu }}^{2}$. For NGC 300 ULX-1, n is −0.8(4). The small magnitude of n indicates that we are likely seeing the secular spin up of the system and the $\ddot{\nu }$ is not strongly contaminated by timing noise or changing Doppler shift from an orbit. An analysis of the long-term spin up in terms of accretion torque theory is presented in Vasilopoulos et al. (2018b).The fitting procedure also yields a pulsed count rate for each observation; these are shown in Figure 3.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

To make more precise frequency measurements, we repeated the likelihood analysis above, with expanded observation segments with a maximum span of 15,000 s (limiting the duration to three NICER orbits). In this stage, there are potential pulse counting ambiguities due to the long gaps between orbits in some observations. However, in nearly all cases, measurements that suffered from an off-by-one counting error between observation subsets were very clearly identified from the short-segment likelihood frequency estimates and the known position of the sidelobes from the observation window function.

These precise frequency measurements, excluding the single observation in 2018 February, are shown in Figure 4, where a constant frequency derivative has been subtracted to allow the fine structure to be seen. This reveals relatively smooth frequency evolution that shows some evidence for several negative frequency steps, which we identify as candidate glitches. However, the sampling and large uncertainties in some of the measurements, coupled with the rapid underlying spin evolution, makes this result tentative, requiring confirmation from the phase-coherent analysis presented in Section 4.

Zoom In Zoom Out Reset image size

We attempted to fit the frequency evolution with a secular evolution up through $\ddot{\nu }$ and three glitches. The residuals to this fit are shown in Figure 5. There are clearly some unmodeled residuals, particularly in the high-precision measurements between glitches 1 and 2, and glitches 2 and 3. This could potentially be modeled as frequency derivative steps at the glitch, or exponentially decaying glitch "recoveries." Both of these phenomena are observed in rotation-powered pulsar glitches, but it is not clear which applies here, so we keep the model as simple as possible. The model parameters are shown in Table 1. We note that these model parameters yield a braking index n = −0.9(2).

Zoom In Zoom Out Reset image size

Table 1.  Incoherent Timing Parameters

Parameter	Value
Epoch (MJD)	58243.515
ν (mHz)	53.3533(7)
$\dot{\nu }$ (Hz s−1)	4.40(1) × 10−10
$\ddot{\nu }$ (Hz s−2)	−3.4(9) × 10−18
${\nu }^{(3)}$ (Hz s−3)	−1.1(3) × 10−24
Glitch 1 Epoch (MJD)	58243.515
Glitch 1 Δν (Hz)	−2.1(1) × 10−5
Glitch 2 Epoch (MJD)	58265.5
Glitch 2 Δν (Hz)	−3.8(2) × 10−5
Glitch 3 Epoch (MJD)	58333.5
Glitch 3 Δν (Hz)	−3.7(5) × 10−5

Download table as:  ASCIITypeset image

3.1. Search for Profile Variations

In most radio pulsar glitches the pulse profile does not change at the epoch of the glitch, though significant profile changes are observed at the epoch of some magnetar glitches (e.g., see Archibald et al. 2013) and radiative changes were seen after a glitch in one high magnetic field pulsar PSR J1119−6127 (Weltevrede et al. 2011). To help constrain the models for the NGC 300 ULX-1 glitches, we looked for profile variations. To do this in a timing-model-independent way, we considered observations where the pulsations were strongly detected (6σ or more). For each segment, we fitted a single Gaussian template to determine phase 0 and generated a set of phases using that as a reference. We then collected all of these phases from six different sections of the NICER monitoring observation, splitting the data at the glitches, and at the first significant data gap after each glitch. We binned these phases into 40-bin histograms, subtracted the mean and divided by the standard deviation, to account for the changing background and evolving pulsed count rate. These six pulse profiles are displayed in Figure 6. We also fitted a three-Gaussian template model to the complete data set and displayed the residuals of each profile to that mean shape. We saw no significant pulse profile changes during our observations.

Zoom In Zoom Out Reset image size

Using the observed frequency evolution, we computed pulse times of arrival (TOAs) for each 150 s segment of NICER data (300 s for segment F, as defined below, where the count rate was significantly lower), using the unbinned maximum likelihood technique described in Ray et al. (2011). The reference phase for each TOA is the peak of the Gaussian pulse. Typical TOA uncertainties are 0.5 s, or about 3% of a pulse period.

Because of the significant red noise, the rapid spin evolution, and the gaps between observations, generating a phase-connected timing model for the full data set is not possible. However, we were able to make phase-connected timing models for several contiguous sections of data (which we label as segments A–F), including spans that cover the first (segment A) and second (segment B) glitches. For later segments (C–F) of observations, we fitted simple models with ν and $\dot{\nu }$. We use Tempo2 (Hobbs et al. 2006) to fit timing models to the TOAs from each segment. The timing model parameters for each segment are presented in Tables 2–4, and the modeled spin-frequency evolution is displayed in Figure 7. The TOA residuals for each model are presented in Figures 8–10.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 2.  Parameters for Glitch 1 (Segment A)

Fit and Data Set
MJD range	58238.9–58249.9
Data span (d)	10
Number of TOAs	385
rms timing residual (s)	0.607
Reduced χ2 value	1.6
Fixed Quantities
R.A., α (hh:mm:ss, J2000)	00:55:04.86
decl., δ (dd:mm:ss, J2000)	−37:41:43.7
Frequency second derivative, $\ddot{\nu }$ (Hz s−2)	0
Epoch of frequency determination (MJD)	58243.5
Glitch epoch (MJD)	58243.509876
Measured Quantities
Pulse frequency, ν (mHz)	53.35285(19)
Frequency first derivative, $\dot{\nu }$ (Hz s−1)	4.358(13) × 10−10
Glitch ${\rm{\Delta }}\nu$	−2.340(20) × 10−5
Glitch ${\rm{\Delta }}\dot{\nu }$	1.41(13) × 10−11
Derived Quantities
Δν/ν	−4.39(4) × 10−4
${\rm{\Delta }}\dot{\nu }/\dot{\nu }$	3.2(3) × 10−2

Download table as:  ASCIITypeset image

Table 3.  Parameters for Glitch 2 (Segment B)

Fit and Data Set
MJD range	58252.5–58278.7
Data span (d)	26
Number of TOAs	235
rms timing residual (μ s)	1.36
Reduced χ2 value	7.0
Fixed Quantities
R.A., α (hh:mm:ss, J2000)	00:55:04.86
decl., δ (dd:mm:ss, J2000)	−37:41:43.7
Frequency second derivative, $\ddot{\nu }$ (Hz s−2)	0
Epoch of frequency determination (MJD)	58265
Glitch epoch (MJD)	58265.25263
Measured Quantities
Pulse frequency, ν (mHz)	54.12969(9)
First derivative of pulse frequency, $\dot{\nu }$ (s−2)	4.1603(17) × 10−10
Glitch Δν	−2.990(14) × 10−5
Glitch ${\rm{\Delta }}\dot{\nu }$	1.24(3) × 10−11
Derived Quantities
Δν/ν	−5.52(2) × 10−4
${\rm{\Delta }}\dot{\nu }/\dot{\nu }$	2.98(6) × 10−2

Download table as:  ASCIITypeset image

Table 4.  Coherent Timing Model Parameters for Segments C–F

Fit and Data Set
Sequence	C	D	E	F
MJD range	58292.6–58296.0	58310.4–58318.4	58327.2–58332.3	58334.2–58353.2
Data span (d)	3.4	8	5.1	19
Number of TOAs	28	89	39	76
rms resid. (s)	0.93	0.88	0.73	2.0
Reduced χ2	3.0	2.3	1.5	15.5
Fixed Quantities
R.A., α	00:55:04.86
decl., δ	−37:41:43.7
$\ddot{\nu }$ (s−3)	0
Epoch (MJD)	58298	58298	58333	58333
Measured Quantities
ν (mHz)	55.3280(8)	55.3295(5)	56.5046(3)	56.46146(13)
$\dot{\nu }$ (Hz s−1)	4.25(3) × 10−10	3.877(3) × 10−10	3.838(9) × 10−10	3.8299(14) × 10−10
Derived Quantities
Δν (Hz)	(D–C) 1.5(9) × 10−6		(F–E) −4.32(3) × 10−5
Δν/ν	2.7(16) × 10−5		−7.66(6) × 10−4

Download table as:  ASCIITypeset image

We considered the possibility of a glitch between segments C and D, but rejected that hypothesis since the extrapolated fits do not appear to require one. For segments E and F, we were able to get a good measurement of the frequency step at the glitch.

The coherent timing models shown in Figures 8–10 clearly demonstrate that there are sudden pulse frequency changes. However, it is conceivable that such features could be just a manifestation of the red noise in the system. To test this, we simulated red noise residuals. Most accretion-powered pulsars show white torque noise (Bildsten et al. 1997), which gives a pulse phase power spectrum with an index of −4, while some are even steeper. We made multiple realizations of residuals with red noise spectral indices of −3, −4, and −5, with normalizations set to match the rms of the observed residuals. The simulations with index −3 sometimes showed cusp-like features, but never with the large difference in scale between the red noise wiggles and the glitch residuals that we see in Figures 8 and 10. For more realistic red noise spectral indices of −4 and −5, the red noise is much smoother and extremely unlikely to be mistaken for a glitch. We therefore conclude that the timing features we observe are best characterized as glitches.

Spectral analysis was performed with XSPEC version 12.10 (Arnaud 1996), and with the NICER response files version 0.06. For these analyses, the data filtering is more stringent than that of the timing analyses presented above. This is in an effort to minimize the non-astrophysical background that affects the spectra. First, we limit the energy range to 0.4–10 keV because of significant optical loading affecting channels below 0.4 keV, especially when NICER is observing at low Sun angles. In addition to the standard filtering described in Section 2, we also exclude the detectors with DET_ID 14, 34, and 54, which often suffer from optical loading. Finally, we exclude periods with the conditions: COR_SAX < 1.5, FPM_OVERONLY_COUNT > 1.0, and ${\mathtt{FPM}}\_{\mathtt{OVERONLY}}\_{\mathtt{COUNT}}\gt \left(1.52\times {\mathtt{COR}}\_{{\mathtt{SAX}}}^{-0.633}\right)$. This last condition is an empirical relation obtained from a large number of "blank sky" pointings to characterize NICER background over a wide range of orbital environments. Together, these remove time intervals with high particle background.

We generated a background spectrum specific to small time intervals using a library of spectra populated with NICER observations of seven blank sky regions, corresponding to the RXTE background regions (Jahoda et al. 2006). The background spectra we obtain are a combination of these blank sky spectra weighted according to space-weather observing conditions common to both the NGC 300 ULX-1 and background-field observations.

To each ObsID spectrum with its corresponding background, we fit a simple absorbed power law to study the spectral evolution over the several weeks of NICER observations of NGC 300 ULX-1. No obvious spectral changes can be determined at the epochs of the glitches. However, large variations of the background count rates estimated for each ObsID may hamper the precise study of the evolution.

Instead, we focus our study on the spectral evolution of the pulsed component only. For each ObsID, we extract spectra from the on-pulse phase range (0.75–1.25) and off-pulse phase range (0.25–0.75; see Figure 6). The off-pulse spectra, used as backgrounds for the on-pulse spectra, therefore permit subtraction of the non-astrophysical background and of the other X-ray sources within the NICER field of view. The typical flux of the pulsed component is f_X ∼ 2 × 10⁻¹² erg s⁻¹ cm⁻² (0.4–10 keV).

We fit each of these on-minus-off spectra with a simple absorbed power law. The best-fit parameters are shown in Figure 11, together with the reduced ${\chi }^{2}$ values. No significant spectral change is observed at the time of glitches 1 and 2. The observations around glitch 3 were shorter and provided low signal-to-noise spectra for a similar study. We therefore conclude that no detectable radiative changes were associated with the glitches.

Zoom In Zoom Out Reset image size

6.1. Characterization of the Glitches

6.2. Nature of the Glitches

The glitches measured in NGC 300 ULX-1 are unique among glitches observed to date. The lack of strong evidence for any radiative change associated with the glitch in NGC 300 ULX-1 suggests a glitch mechanism internal to the neutron star, rather than a magnetospheric mechanism. Even in the case of magnetar glitches, where only some glitches are accompanied by radiative changes, they are thought to be caused by a source within the neutron star (Dib & Kaspi 2014). The internal mechanism that produces glitches in nonaccreting pulsars is thought to be due to a superfluid in the pulsar's inner crust and possibly core. The free neutrons that are present there are expected to be in a superfluid state because the typical crust temperature is well below the critical temperature for superfluidity. Unlike normal matter, superfluid matter rotates by forming vortices whose areal density determines the spin rate of the superfluid. These vortices are normally pinned to the lattice of crust nuclei, and as a result, while the rest of the star slows down owing to electromagnetic dipole radiation, the superfluid does not. A lag develops between the stellar spin rate and that of the superfluid. When this lag reaches a critical value, superfluid vortices unpin and transfer angular momentum to the rest of the star, causing the stellar spin rate to increase, i.e., a spin-up glitch (Anderson & Itoh 1975; see, e.g., Haskell & Melatos 2015, for a review; see left panel of Figure 12).

Zoom In Zoom Out Reset image size

The above mechanism can still apply to the case of a spin-down glitch (see, e.g., Pines et al. 1980; see right panel of Figure 12), such as that seen here for NGC 300 ULX-1. Here, the pulsar spin frequency increases with time due to accretion (pre-glitch $\dot{\nu }=4\times {10}^{-10}\,\mathrm{Hz}\,{{\rm{s}}}^{-1}$), with a characteristic spin-up timescale ${\tau }_{{\rm{c}}}\equiv \nu /2\dot{\nu }=2\,\mathrm{yr}$. Because of this short timescale, we might expect a lag to build up quickly, where in our case most of the star spins faster than the superfluid. Once the lag is too great, the superfluid receives angular momentum from the rest of the star, and the star undergoes a spin-down glitch. This can be shown from the calculation of e.g., Andersson et al. (2012), by simply replacing the torque due to dipole radiation with that due to accretion. Such an interpretation provides a constraint on the superfluid fraction in the neutron star (Link et al. 1999; Andersson et al. 2012; Chamel 2013) and even potentially a determination of the neutron star mass (Ho et al. 2015, 2017). The former yields ${I}_{\mathrm{sf}}/I\approx ({m}_{{\rm{n}}}^{\mathrm{eff}}/{m}_{{\rm{n}}})$ $(-{{\rm{\Sigma }}}_{i}{\rm{\Delta }}{\nu }^{i}/\nu )$ $(2{\tau }_{{\rm{c}}}/{t}_{\mathrm{obs}})$ ≈ $0.02({m}_{{\rm{n}}}^{\mathrm{eff}}/{m}_{{\rm{n}}})$ $(116\,{\rm{d}}/{t}_{\mathrm{obs}})$, where the superfluid and total moments of inertia are I_sf and I, respectively, the average effective neutron mass due to entrainment is m^eff_n/m_n (∼4.3; Chamel 2012), the summation is over each glitch, and the observation time is t_obs. The large glitch amplitudes ($| {\rm{\Delta }}\nu | \approx 20-40$ μHz) and frequency derivative ($\dot{\nu }=4\times {10}^{-10}\,{{\rm{Hz}}{\rm{s}}}^{-1}$) of NGC 300 ULX-1 are comparable to those of the frequent large glitching pulsars Vela and PSR J0537−6910 and ${\rm{\Delta }}\dot{\nu }\approx {10}^{-11}\,{{\rm{Hz}}{\rm{s}}}^{-1}$ are among the largest measured (Espinoza et al. 2011; Antonopoulou et al. 2018; Ferdman et al. 2018). This could suggest that NGC 300 ULX-1 has a similarly high glitch activity (Espinoza et al. 2011; Yu et al. 2013; Fuentes et al. 2017) if the underlying glitch mechanism is the same.

The rise time of a glitch is uncertain, due to unknown properties of superfluidity, but is thought to be >8 P (Sidery et al. 2010; Sourie et al. 2017). Such a limit is in accord with glitches seen in the Vela pulsar (Palfreyman et al. 2018) and suggests a short rise time of >150 s for NGC 300 ULX-1, well within the ∼12 hr of our nearest post-glitch TOA (see Section 3).

6.3. Glitch Recurrence Timescale

6.4. No Evidence for an Orbit

6.5. Magnetic Field and Energetics

In this article, we presented evidence for two spin-down glitches (and a tentative third one) during the spin up of the ULX source NGC 300 ULX-1. Our NICER monitoring campaign allowed us, using phase-coherent timing analyses, to determine with high confidence that the magnitudes of these two glitches are the largest among known radio pulsar glitches. While the third glitch had an even larger magnitude (Δν = −43 μHz), its identification is somewhat less secure, as gaps in monitoring precluded precise determination of the glitch epoch and amplitude. We argue that these spin-down glitches are similar in nature to the spin-up glitches observed in spinning-down radio pulsars, i.e., the glitches are caused by sudden angular momentum transfer between a superfluid component and the rest of the neutron star. This is supported by the lack of detected radiative changes, which indicates a process internal to the neutron star. Finally, despite our dense monitoring campaign, our timing solution did not permit determination of the orbit of the binary system, either because the system's inclination is near 0 ° or because the Doppler amplitude of a long period orbit (hundreds of days) is hidden in the frequency evolution ($\dot{\nu }$ and $\ddot{\nu }$) of the accretion torques in this system. Continued dense monitoring would permit studying the spin period evolution and accretion mechanisms of this intriguing source, as well as detecting more spin-down glitches to characterize the superfluid components inside the neutron star.

We thank the anonymous referee for useful comments and suggestions that improved this article. We thank Nishanth Anand (Thomas Jefferson High School for Science and Technology) for his work on the analysis for this paper during his SEAP summer internship at NRL. We thank Faith Huynh (MIT) for developing excellent NICER quicklook plotting software used in this work. We also thank Breanna Binder for helpful discussions. This work was supported by NASA through the NICER mission and the Astrophysics Explorers Program. NICER work at NRL is supported by NASA. S.G. is supported by the Centre National d'Etudes Spatiales (CNES). W.C.G.H. is partially supported by the Science and Technology Facilities Council (STFC) in the UK. This research has made use of data and software provided by the High Energy Astrophysics Science Archive Research Center (HEASARC), which is a service of the Astrophysics Science Division at NASA/GSFC and the High Energy Astrophysics Division of the Smithsonian Astrophysical Observatory. This research has also made use of the NASA Astrophysics Data System (ADS) and the arXiv.

Facility: NICER. -

Software: astropy (Astropy Collaboration et al. 2013), PINT (https://github.com/nanograv/pint), HEAsoft (https://heasarc.nasa.gov/lheasoft/), Tempo2 (Hobbs et al. 2006), PRESTO (https://www.cv.nrao.edu/~sransom/presto/).