Issue 
A&A
Volume 641, September 2020



Article Number  L4  
Number of page(s)  11  
Section  Letters to the Editor  
DOI  https://doi.org/10.1051/00046361/202038378  
Published online  01 September 2020 
Letter to the Editor
Extreme intrahour variability of the radio source J1402+5347 discovered with Apertif^{⋆}
^{1}
ASTRON, the Netherlands Institute for Radio Astronomy, Oude Hoogeveensedijk 4, 7991 PD Dwingeloo, The Netherlands
email: oosterloo@astron.nl
^{2}
Kapteyn Astronomical Institute, University of Groningen, Postbus 800, 9700 AV Groningen, The Netherlands
^{3}
Astro Space Center of Lebedev Physical Institute, Profsoyuznaya Str. 84/32, 117997 Moscow, Russia
^{4}
Astronomisches Institut der RuhrUniversität Bochum (AIRUB), Universitätsstrasse 150, 44780 Bochum, Germany
^{5}
Department of Astronomy, University of Cape Town, Private Bag X3, Rondebosch 7701, South Africa
^{6}
Department of Physics, Virginia Polytechnic Institute and State University, 50 West Campus Drive, Blacksburg, VA 24061, USA
^{7}
Anton Pannekoek Institute, University of Amsterdam, Postbus 94249, 1090 GE Amsterdam, The Netherlands
^{8}
CSIRO Astronomy and Space Science, Australia Telescope National Facility, PO Box 76, Epping, NSW 1710, Australia
^{9}
Sydney Institute for Astronomy, School of Physics, University of Sydney, Sydney, NSW 2006, Australia
^{10}
Tricas Industrial Design & Engineering, Hanzelaan 95b, 8017 JE Zwolle, The Netherlands
^{11}
Center for Information Technology, University of Groningen, Postbus 11044, 9700 CA Groningen, The Netherlands
Received:
8
May
2020
Accepted:
17
August
2020
The propagation of radio waves from distant compact radio sources through turbulent interstellar plasma in our Galaxy causes these sources to twinkle, a phenomenon called interstellar scintillation. Such scintillations are a unique probe of the microarcsecond structure of radio sources as well as of the subAUscale structure of the Galactic interstellar medium. Weak scintillations (i.e. an intensity modulation of a few percent) on timescales of a few days or longer are commonly seen at centimetre wavelengths and are thought to result from the lineofsight integrated turbulence in the interstellar plasma of the Milky Way. So far, only three sources were known that show more extreme variations, with modulations at the level of some dozen percent on timescales shorter than an hour. This requires propagation through nearby (d ≲ 10 pc) anomalously dense (n_{e} ∼ 10^{2} cm^{−3}) plasma clouds. Here we report the discovery with Apertif of a source (J1402+5347) showing extreme (∼50%) and rapid variations on a timescale of just 6.5 min in the decimetre band (1.4 GHz). The spatial scintillation pattern is highly anisotropic, with a semiminor axis of about 20 000 km. The canonical theory of refractive scintillation constrains the scattering plasma to be within the Oort cloud. The sightline to J1402+5347, however, passes unusually close to the B3 star Alkaid (η UMa) at a distance of 32 pc. If the scintillations are associated with Alkaid, then the angular size of J1402+5347 along the minor axis of the scintels must be smaller than ≈10 μas, yielding an apparent brightness temperature for an isotropic source of ≳10^{14} K.
Key words: scattering / ISM: clouds
Light curves are only available at the CDS via anonymous ftp to cdsarc.ustrasbg.fr (130.79.128.5) or via http://cdsarc.ustrasbg.fr/vizbin/cat/J/A+A/641/L4
© ESO 2020
1. Introduction
Intrahour variability (IHV) refers to quasirandom intensity modulation that is similar to canonical interstellar scintillation, but on vastly shorter timescales. It is an extreme manifestation of radiowave propagation through turbulent plasma, but the nature of the plasma clouds causing the fast variability has remained unclear so far. It is generally thought to be caused by nearby (d ∼ 10 pc) turbulent plasma (Bignall et al. 2007; de Bruyn & Macquart 2015). Because the scattering strength of plasma is weaker if the plasma clouds are at small distances, the plasma densities inferred from IHV observations are orders of magnitude higher than that in the ambient ionised interstellar gas. The formation and survival of such overpressured clouds is at the heart of the mystery. No multiwavelength counterparts to the plasma have yet been found.
A recent development is the statistical association of three known IHVs with nearby hot stars by Walker et al. (2017) suggested that the scattering is caused by the ionised sheath of tiny selfgravitating molecular clumps similar to the cometary knots seen in the Helix nebula, for example. If true, this would be extremely interesting as these clouds would constitute a new population of objects within the interstellar medium. Alternatively, Pen & King (2012) and Pen & Levin (2014) have refined the model of Goldreich & Sridhar (2006), which suggests that the extreme scattering phenomena are caused by thin plasma sheets (which have a basis in magnetohydrodynamics, MHD, theory) oriented along the sightline. The implication here being that there is nothing inherently unusual about the intervening plasma save a fortuitous orientation.
Only very few IHVs are known. More such objects need to be discovered to better understand their statistical properties and for more secure empirical associations. More cases with extreme properties need to be studied to test the limits of proposed theories. The recent spate of widefield survey radio telescopes operating in the GHz band, such as the Australian Square Kilometre Array Pathfinder (ASKAP) and the Aperture Tile In Focus (Apertif) project on the Westerbork Synthesis Radio Telescope (WSRT), are in an excellent position to provide larger samples. Here we report on the first IHV discovered by the recently started survey with the Apertif phasedarray frontends on the WSRT.
2. Discovery and followup
We serendipitously discovered large and rapid scintillations in the radio source J1402+5347 (RA_{J2000} 14^{h} 02^{m} , Dec_{J2000} 53° 47′ 11″) with the Apertif phasedarray feed system on the WSRT in commissioning data taken on April 9, 2019. Because of the eastwest orientation of the WSRT, highly variable sources show prominent artefacts in the form of radial spikes in longexposure images, which makes them easy to identify (Fig. 1). The light curve extracted for the source (Fig. 2; details in Appendix A) shows strong fluxdensity variations that are about 38% of the mean value. Archival WSRT observations from 2013 of the same field at the same frequency showed no significant variability in J1402+5347. No counterpart to J1402+5347 is visible in optical images, but an object is present at the location of J1402+5347 in images produced by the WISE satellite (Cutri et al. 2012). The lack of an optical counterpart together with the colours of this object in the WISE data suggests that the background object is a dustobscured active galactic nucleus (AGN) at fairly high redshift (z ≳ 1; Jarrett et al. 2017).
Fig. 1. Discovery image of J1402+5347 (marked by the blue lines) made from the observation on April 9, 2019. The nearby spiral galaxy M 101 lies to the north. J1402+5347 was identified as a highly variable source because of the prominent radial spikes that are centred on it. 
Fig. 2. Light curves of J1402+5347 observed over a year with a monthly cadence. The observation dates are indicated. The annual modulation in the variation rate is apparent. We observe two standstills in August 2019 and November 2019, while the most rapid variations occurred in April 2019 and March 2020. 
Followup observations taken in May 2019 also showed strong variability at about 30%. Further observations in June and July 2019 showed weaker and slower variability (Fig. 2). An annual modulation of the scintillation properties has been observed in the other known IHVs (DennettThorpe & de Bruyn 2003; Walker et al. 2009; Bignall et al. 2019) and is caused by the annual change in the relative velocity of the Earth with respect to the plasma clouds that is due to the orbital motion of the Earth around the Sun. Motivated by this, we continued monitoring the source with a cadence of about one month. The intensity and rate of the variations continued to decrease and J1402+5347 virtually ceased to vary in August. Following a brief period of slow variability in September and October, another standstill was observed in November. After November, the rate and intensity of the variations increased steadily, and by March 2020, they even exceeded the high levels seen in the discovery observation of April 2019.
3. Timescale analysis
The scintillation geometry and physical properties of the plasma clouds causing the scintillation can be determined from the annual variation in the scintillation timescale. We computed the scintillation timescale for each observation using two methods that are widely used in the analysis of timeseries. The first method determines the decorrelation timescale from the noisebiascorrected ACF directly computed from the light curves (Appendix B). The temporal lag at which the ACF of the variations falls below a specified threshold is taken as the scintillation timescale. For easy comparison with previous work (DennettThorpe & de Bruyn 2003; Bignall et al. 2019), we define the threshold as 1/e times the peak value. This method is more accurate as it makes minimal assumptions. However, it is only effective for epochs of rapid variability as it requires several scintels to be observed within a single observation run (DennettThorpe & de Bruyn 2003). For epochs near the standstills, we used a method that is commonly called Gaussian process (GP) regression (Appendix B). This method models the observed variations as a GP with a judiciously chosen functional form for the temporal covariance function. Although this method assumes an analytic form of the covariance function, it can effectively determine the decorrelation timescale and uncertainties even in epochs near the standstills (Bignall et al. 2019). The decorrelation timescales from both methods are tabulated in Table 1 and agree with each other within the errors for all the epochs to which both methods could be applied.
Annual modulation statistics.
The scintillation timescale of just 6.5 min in March 2020 and the similarly short timescale in the April 2019 observations represent the fastest interstellar scintillation reported to date. This is remarkable because, compared to previous studies, our observations are at lower frequency for which normally the scintillations are significantly slower.
This timescale directly relates the physical size of the scintels to the transverse speed between the Earth and the scattering plasma. A singleepoch measurement cannot simultaneously constrain the scintel size and the relative speed. However, the annual modulation in the scintillation timescale can be used to estimate both parameters. In particular, the two standstills critically constrain both the screen velocity and the geometry of the scintels. One annual standstill can be obtained when the screen velocity on the ecliptic plane equals the Earth’s orbital velocity. The second standstill requires the scintels to be anisotropic such that the standstill occurs when the relative velocity between the Earth and the plasma cloud is along the major axis of the scintels. We fit a simple model with onedimensional scintels to the observed annual modulation, as suggested by Walker et al. (2009, see Appendix C). The annual variation in this model is fully specified by three parameters: (a) the semiminor axis of the scintels, a_{⊥}, (b) the systemic velocity of the scattering plasma along the minor axis of the scintels, v_{⊥}, and (c) the orientation of the major axis on the plane of the sky, θ_{R} (measured from north to east). Figures 3 and 4 show the posterior distributions of the model parameters and the bestfitting annual modulation curve, respectively. We constrain the semiminor axis of the scintels to be just a_{⊥} = 2.17(0.08)×10^{4} km. We note that the reduced chisquared for the bestfitting model is close to unity, which implies that the onedimensional scintel model sufficiently and parsimoniously captures the annual variation in the data. To understand the degree of anisotropy implied by the data, we also computed the posterior distribution for a twodimensional scintillation model (Appendix C). The twodimensional model places a lower limit on the ratio of the long and short axes of the scintels, R, of log_{10}(R) > 2.12. As expected from the high anisotropy, the data do not constrain the space velocity of the scattering screen parallel to the long axis of the scintels.
Fig. 3. Posterior distribution of the model parameters describing the annual modulation of the scintillation rate. The panels on the diagonal are the marginalised onedimensional posterior distributions. The dashed black lines are drawn at 0.16, 0.5, and 0.84 quartiles. The annotated text gives the bestfit parameters and ±1σ bounds implied by the posterior distributions. 
Fig. 4. Scintillation rate determined through a Gaussian regression analysis (blue points with ±1σ errors); the maximum likelihood model (black curve) is overplotted. The annotated text gives the model parameters and the reduced χ^{2} of the fit. The phase at the vernal equinox is defined to be zero, with phase increasing with time. Magenta crosses show the inverse of the scintillation timescale measured directly from the autocorrelation function (ACF) of the light curves during epochs of sufficiently rapid variations. 
4. Discussion and outlook
Incoherent synchrotron sources are generally expected to only display refractive interstellar scintillation. Refractive scintels of a pointlike source will have a size equal to or larger than that of the first Fresnel zone: , where λ is the wavelength and d is the distance to the scattering screen (Goodman & Narayan 2006). If the measured semiminor axis of the scintels is associated with the Fresnel scale, the distance to the scattering plasma is d = 0.24 pc, which is within the sphere of gravitational influence of the Sun. This would place the scattering plasma within the canonically defined Oort cloud.
If J1402+5347 is unusually compact, then it may display diffractive interstellar scintillations, as is usually seen in pulsars. In this case, the scintel size is associated with the diffractive scale. Because the diffractive scale can be much smaller than the Fresnel scale, the data admit significantly larger screen distances in the diffractive regime. Regardless, the length scale of scintels cannot be smaller than the projected size of the source at the scattering screen. We consider for simplicity a circular source of angular radius θ_{src}. The requirement d θ_{src} < a_{⊥} places a constraint on the apparent brightness temperature of the source of T_{b} > 10^{11}(d/pc)^{2} K.
Walker et al. (2017) have shown a statistical astrometric association between sightlines to sources showing intrahour variability and hot stars such as Vega, Spica, and Alhakim. They also reported that the long axis of the scintels of IHVs preferentially points towards the associated star. We searched the HIPPARCOS catalogue to verify if the sightline to J1402+5347 passes anomalously close to any known O, B, or Atype star (see Fig. 5). The B3 star Alkaid (η UMa, HIP 67301) at a distance of d = 32 pc is unusually close to the sightline towards J1402+5347. The volume density of Btype stars in the HIPPARCOS catalogue that lie within a 50degreewide cone around this sightline is about 4.3 × 10^{−5} pc^{−3}. The impact parameter of the sightline with Alkaid is about 2.7 pc and the angle between the scintel long axis and the starward direction for a pairing of J1402+5347 and Alkaid is about 9 ± 3 deg. The probability of finding a Btype star by chance with the observed impact parameter or closer and with the observed angular alignment or better is about 0.2%. In addition, the space velocity components of Alkaid with respect to the Solar System barycentre are −17.8 km s^{−1} and −2.2 km s^{−1} along the RA and Dec axes, respectively. This is about 9 km s^{−1} different from the space velocity component of the scattering plasma as constrained by our observations, although, as pointed out by Walker et al. (2017), a clear agreement between space velocities may be precluded because of the free expansion of the scattering plasma, and to a lesser extent, because of projection effects. Regardless, an association of the scattering plasma with Alkaid sets d = 32 pc and places a rather extreme lower limit on the apparent brightness temperature of J1402+5347 of T_{b} > 10^{14} K. This limit can be slightly ameliorated if we postulate that the radio source itself is extended precisely along the major axis of the scintels.
Fig. 5. Left panel: line of sight and transverse distance of stars in the HIPPARCOS catalogue with respect to the sightline towards J1402+5347. The stars are colourcoded by their spectral type (legend in the annotated text). Right panel: same as the left panel, but the stars are plotted in the plane of the sky. The central black square shows J1402+5347. The solid and dashed blue lines mark the orientation of the major and minor axes of the scintels, respectively. 
An association with Alkaid also has critical implications for the variance in the scattering properties of the Galactic plasma at large. The angular broadening expected from the lineofsight integrated Galactic plasma (unrelated to the scattering screen discussed here) for this sightline is estimated to be about 0.36 mas according to the widely used NE2001 model (Cordes & Lazio 2002). The brightness temperature necessary for an association with Alkaid, however, requires that the source diameter does not exceed about 8.6 μas, which is 40 times smaller than the anticipated angular broadening. An association with Alkaid therefore suggests that the bulk of the scattering in the Galactic interstellar medium must occur in highly localised clumps so that this strong sightlinetosightline variation can be accounted for.
Whether the scattering plasma is located in the Oort cloud or is associated with Alkaid may be determined by extending theoretical analyses, such as the analysis of Goodman & Narayan (2006), to the case of highly anisotropic scintels, and by broadband observations of J1402+5347 at higher frequencies. If the observed scattering at 1.4 GHz is indeed diffractive, two effects must be seen in such observations: (a) the scintillations will decorrelate on frequency scales Δν/ν ≪ 1, and (b) the transition to weak scattering must occur at frequencies that are several times higher than 1.4 GHz, and similarly strong lightcurve variations must then extend to these high frequencies. If the scattering is associated with the Oort cloud, then the observed variations are due to refractive scintillations in the weakscattering regime. We therefore expect to observe broadband variations (Δν/ν ∼ 1) at higher frequencies with a continuously declining level of fractional variations with increasing frequency. The observations of April 9, 2019, show that the variations at the upper and lower end of the Δν/ν ≈ 15% observing band are very similar, with a Pearson correlation coefficient of 0.84 (Appendix D), suggesting refractive variations and hence that the scattering plasma is located in the Oort cloud. However, further theoretical investigation on the spectral decorrelation properties of highly anisotropic scintels are required to clearly determine the distance to the scattering screen.
We end by noting that because of the larger apparent source sizes at lower frequencies, fast intrahour variability at decimetre wavelengths will be significantly biased towards intervening plasma clouds that are nearby. The discovery of rapid scintillation in J1402+5347 in early Apertif data bodes well for ongoing widefield surveys in the decimetre band. Based on the discovery of J1402+5347 and the sky area covered by the Apertif survey at the time, we expect to discover approximately ten more such systems, and indeed, several new IHV sources have more recently been identified in Apertif data (Oosterloo et al., in prep.).
Acknowledgments
We dedicate this Letter to the memory of Ger de Bruyn and of J.P. Macquart. HKV thanks Mark Walker for discussions. Credit is also due to Mark Walker for pointing out the possible association with Alkaid. We thank the anonymous referee for a meticulous reading of the manuscript which helped us to improve our paper, and for rectifying some errors in notation. EAKA is supported by the WISE research programme, which is financed by the Netherlands Organisation for Scientific Research (NWO). YM acknowledges funding from the European Research Council under the European Union’s Seventh Framework Programme (FP/20072013)/ERC Grant Agreement No. 617199. LCO acknowledges funding from the European Research Council under the European Union’s Seventh Framework Programme (FP/20072013)/ERC Grant Agreement No. 617199. JMvdH acknowledges funding from the European Research Council under the European Union’s Seventh Framework Programme (FP/20072013)/ERC Grant Agreement No. 291531 (“HIStoryNU”), JvL acknowledges funding from the European Research Council under the European Union’s Seventh Framework Programme (FP/20072013)/ERC Grant Agreement No. 617199 (“ALERT”), and from Vici research programme “ARGO” with project number 639.043.815, financed by the Netherlands Organisation for Scientific Research (NWO). This work makes use of data from the Apertif system installed at the Westerbork Synthesis Radio Telescope owned by ASTRON. ASTRON, the Netherlands Institute for Radio Astronomy, is an institute of the Dutch Science Organisation (De Nederlandse Organisatie voor Wetenschappelijk Onderzoek, NWO).
References
 Bignall, H. E., Jauncey, D. L., Lovell, J. E. J., et al. 2003, ApJ, 585, 653 [Google Scholar]
 Bignall, H. E., Jauncey, D. L., Lovell, J. E. J., et al. 2007, Astron. Astrophys. Trans., 26, 567 [Google Scholar]
 Bignall, H., Reynolds, C., Stevens, J., et al. 2019, MNRAS, 487, 4372 [Google Scholar]
 Cordes, J. M., & Lazio, T. J. W. 2002, ArXiv eprints [arXiv:astroph/0207156] [Google Scholar]
 Cutri, R. M., Wright, E. L., Conrow, T., et al. 2012, Explanatory Supplement to the WISE AllSky Data Release Products [Google Scholar]
 de Bruyn, A. G., & Macquart, J. P. 2015, A&A, 574, A125 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 DennettThorpe, J., & de Bruyn, A. G. 2003, A&A, 404, 113 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 ForemanMackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306 [CrossRef] [Google Scholar]
 ForemanMackey, D., Agol, E., Ambikasaran, S., & Angus, R. 2017, AJ, 154, 220 [NASA ADS] [CrossRef] [Google Scholar]
 Goldreich, P., & Sridhar, S. 2006, ApJ, 640, L159 [Google Scholar]
 Goodman, J., & Narayan, R. 2006, ApJ, 636, 510 [Google Scholar]
 Jarrett, T. H., Cluver, M. E., Magoulas, C., et al. 2017, ApJ, 836, 182 [Google Scholar]
 Pen, U.L., & King, L. 2012, MNRAS, 421, L132 [Google Scholar]
 Pen, U.L., & Levin, Y. 2014, MNRAS, 442, 3338 [Google Scholar]
 Walker, M. A., de Bruyn, A. G., & Bignall, H. E. 2009, MNRAS, 397, 447 [Google Scholar]
 Walker, M. A., Tuntsov, A. V., Bignall, H., et al. 2017, ApJ, 843, 15 [Google Scholar]
Appendix A: Interferometric data reduction
The observations of J1402+5347 presented here were performed with the Apertif phasedarray feed system that was recently installed on the WSRT. With Apertif, 40 partially overlapping beams, each with a full width at half maximum of about 35 arcmin, are formed on the sky. The data from each beam are used independently to image the region covered by each beam using standard aperture synthesis. The total area imaged by combining the images from all beams is about with a spatial resolution of 12″ × 17″ at the declination of J1402+5347. The observing band was 122 MHz wide, divided into 156 channels, centred on 1.365 GHz. J1402+5347 was detected in three overlapping beams. Each observation had a duration of 11 to 11.5 h. The calibration of the data from each beam followed standard procedures, involving excision of radiofrequency interference, flux, bandpass, and crosscalibration using 3C147, followed by selfcalibration. The resulting images had small directiondependent imaging errors, which were corrected for using standard peeling techniques. This applies in particular to the strong source 3C295, located SW of J1402+5347: because of the large distance from the beam centres, it showed large directiondependent effects. The resulting images are noise limited with a noise level of about 30 μJy beam^{−1} that varies slightly between observations.
To construct the light curves, the source model resulting from the selfcalibration was subtracted from the visibilities, except for the source model for J1402+5347. This created a visibility data set that contained only emission from J1402+5347. The phase centre of this visibility data set was shifted to the position of J1402+5347 so that the phases of all visibilities were zero. Finally, for each time stamp, the visibilities of all baselines and frequencies were summed to obtain the flux density of J1402+5347 for each time stamp. After primary beam correction, this yields the light curves shown in Fig. 2.
To verify the accuracy of the light curves, the same procedure was followed for J1403+5349, which has a very similar flux density, 9 arcminutes W of J1402+5347. An example of such a comparison is shown in Fig. A.1.
Fig. A.1. Light curves derived for April 9, 2019, of the variable source J1402+5347 and the steady comparison source J1403+5349, 9 arcmin west of J1402+5347, illustrating the accuracy of the method for constructing the light curves. The horizontal axis shows the time since the start of the observation. 
Appendix B: Determination of the scintillation timescale
We computed the scintillation timescale using the autocorrelation function (ACF) and Gaussian process modelling (GP). For the estimates based on the ACF, we followed a procedure that closely mimics that presented by Bignall et al. (2003). We first subtract the mean flux density from the light curves to obtain a time series, x_{i}, with a sampling interval Δt. For each pair ij we compute the product x_{i}x_{j}. The list of products is then binned on a grid of size Δt, and all products falling into the same bin were averaged. Thermal noise only biases the zeroth bin as it is uncorrelated between lightcurve samples. This bias, equal to the variance of the noise, is estimated to be equal to half the variance of the difference between successive samples in the time series. The 1/e timescale, t_{e}, is estimated using a straightline fit to a segment of the ACF containing five samples, centred on the sample that is closest to 1/e. The error on the scintillation timescale is given by , where N is the number of observed cycles of variation on timescale t_{e} and ϵ is a factor of about unity that depends on the exact profile of the ACF around t_{e}. If the total duration of observation is T (∼11.5 h in our case), then N = (T/t_{e}). The factor ϵ was estimated with numerical Monte Carlo techniques (DennettThorpe & de Bruyn 2003; Bignall et al. 2003) to be in the range 0.25 − 0.9. We conservatively assumed ϵ = 1. The timescale and its errors thus computed are only accurate when several scintels are observed, which is not the case during the standstills. An example ACF is plotted in Fig. B.1.
Fig. B.1. Numerical ACF computation: mean subtracted light curve (top panel) and the ACF (bottom panel) for the April 9, 2019, light curve. 
Our second approach was to model the temporal scintillation as a GP with a covariance function given by a damped exponential: K(Δt) = K_{0}exp(cΔt)cos(dΔt), where Δt is the temporal lag and K_{0}, c, and d are model parameters (Bignall et al. 2019). This choice for the covariance function is motivated by the shape of the ACF, which shows a steep rolloff at small lags, follows by oscillatory behaviour due to the dominant spectral mode in the light curves (see Fig. B.3) We assume that the scintillation is a stationary process, which means that the parameters remain the same at different epochs, and the variation in the scintillation timescale can be modelled by a single epochdependent scale factor, k, by substituting Δt → kΔt. This allows us to estimate the scintillation timescale and the associated uncertainty even for epochs close to the standstill by employing our knowledge of the covariance function from epochs of rapid variability. We determined the parameters K_{0}, c, and d from the April 2019 light curve alone instead of using a fit using the entire annual data set (Bignall et al. 2019). With these parameter values fixed, we determined the temporal scale factor k for each epoch with a singleparameter Gaussian regression.
All data were mean subtracted prior to regression analysis, and a noise variance term (estimated from the time differences) was added to the diagonal elements of the covariance matrix. For the April 2019 data set, we fixed K_{0} to equal the total variance in the scintillations: K_{0} = 3.477105 mJy^{2}, and we estimated c and d. The bestfit estimates were log(c) = − 3.80 ± 0.10 and log(d) = − 1.86 ± 0.05. The posterior distributions for these parameters were computed using the emcee software (ForemanMackey et al. 2013) and are plotted in Fig. B.2. The numerically computed ACF and the bestfit kernel are plotted together for comparison in Fig. B.3.
Fig. B.2. Bestfitting kernel: numerically computed ACF overplotted with the bestfit Gaussian kernel for the April 9, 2019, light curve. The decorrelation timescale of the Gaussian kernel is indicated with the dashed black line. The quasisinusoidal oscillations at long delay are not ideally captured by the analytic form of the kernel, However, the decay to the 1/evalue at short delays is well captured and is most relevant to estimating scintillation rates. 
Fig. B.3. Gaussian process kernel parameters: posterior distribution of GP kernel parameters for the April 9, 2019, data set. The total variance K_{0} was fixed to equal the noisebias corrected variance of the light curve, and only the scale parameters c and d were estimated. The dashed lines are placed at the 16, 50, and 84^{th} quantiles of the marginalised onedimensional distributions. The fractional error in the estimation of parameter c is about 10%. 
Next, we held K_{0}, c and d fixed at the bestfit values quoted above and used the celerite (ForemanMackey et al. 2017) package to fit a single scale parameter k. The 1/e crossing timescale for the April 2019 epoch is largely set by the value of parameter c, which is determined with a formal error of about 10%. To account for this uncertainty in our error budget, we added a 10% error in quadrature to the formal errors on the timescale from the GP regression analysis.
Throughout our analysis, the likelihood computations were made using routines from the celerite package (ForemanMackey et al. 2017), the bestfit values were determined using scipy.optimize with the LBFGSB method, and the posterior distributions were sampled using the emcee package (ForemanMackey et al. 2013).
Appendix C: Annual modulation
The timescale of a temporal scintillation caused by anisotropic scattering is given by (Walker et al. 2009)
where R = a_{}/a_{⊥} is the anisotropy of the scintels, a_{⊥} and a_{} are the length scales of scintels along the minor and major axes, respectively, and v_{⊥ rel} and v_{ rel} are the relative velocity between the Earth and the scintels along the minor and major axes of the scintels, respectively. For highly anisotropic scintels (Walker et al. 2009), we can set R → ∞ and obtain t_{s} = a_{⊥}/v_{⊥ rel}. The annual modulation profile can be fit with a threeparameter model. The parameters are (i) the semiminor axis of the scintels, a_{⊥}, (ii) the orientation of the major axis of the scintels in the plane of the sky, θ_{R} ∈ (0, 2π], measured from north towards east, and (iii) the screen velocity in the plane of the sky perpendicular to the major axis of the scintels v_{⊥} ≥ 0 (defined to be along θ_{R} − π/2). We carried out our fits in a tangent plane perpendicular to the sightline with cardinal axes along the directions of increasing RA () and Dec (). If is the velocity of the Earth projected onto this plane, then the relative velocity perpendicular to the major axis of the scintels is given by . We fit this equation to the scintillation timescales tabulated in Table 1 (column “GP scintillation rate”). We carried out the fits to the scintillation rate instead of the scintillation timescale to avoid infinities in the model. First, we carried out a bruteforce search for the maximum likelihood value on a coarse threedimensional parameter grid spanned by 10^{8} cm < a_{⊥} < 5 × 10^{10} cm, 0 km s^{−1} < v_{⊥} < 100 km s^{−1}, and 0 < θ_{R} < 2π. The likelihood peaked at about a_{⊥} ≈ 2 × 10^{9} cm, v_{⊥} ≈ 22 km s^{−1}, and θ_{R} ≈ 0.7. We then used the emcee software (ForemanMackey et al. 2013) to obtain the final fit values and their formal errors. We used 200 randomly initialised walkers, each walking 2000 steps. We imposed a flat prior on the parameters in the range 2 × 10^{8} cm < a_{⊥} < 2 × 10^{10} cm, 15 km s^{−1} < v_{⊥, R} < 35 km s^{−1}, and 0.1 < θ_{R} < 1.1. The posterior parameter distributions were computed and plotted from the emcee samples using the corner package.
We also ran a separate MCMC simulation to sample the posterior distribution for the twodimensional scintillation model given in Eq. (C.1). The additional parameters in the twodimensional model and v_{} and R. We placed the prior bounds v_{} ∈ [0, 100) km s^{−1} and log_{10}R ∈ [0, 4), with the relative velocity parallel to the scintel long axis given by . The posterior distributions are given in Fig. C.1.
Appendix D: Broadband nature of scintels
To investigate a possible frequency dependence of the intensity fluctuations, we constructed separate light curves for the upper and lower halves of the observing band for the observations on April 9, 2019. These separate light curves are shown in Fig. D.1, where the blue line represents the light curve for the frequency range 1311−1408 MHz and the red line shows the range 1408−1505 MHz. It is clear from this figure that the light curves are very similar, showing that there is no frequency structure with the observing band. The Pearson correlation coefficient between the two light curves is 0.84.
Fig. D.1. Separate light curves for the lower and upper halves of the observing band of the discovery observation on April 9, 2019. The blue line is the light curve for the frequency range 1311−1408 MHz, and the red line shows the range 1408−1505 MHz. 
All Tables
All Figures
Fig. 1. Discovery image of J1402+5347 (marked by the blue lines) made from the observation on April 9, 2019. The nearby spiral galaxy M 101 lies to the north. J1402+5347 was identified as a highly variable source because of the prominent radial spikes that are centred on it. 

In the text 
Fig. 2. Light curves of J1402+5347 observed over a year with a monthly cadence. The observation dates are indicated. The annual modulation in the variation rate is apparent. We observe two standstills in August 2019 and November 2019, while the most rapid variations occurred in April 2019 and March 2020. 

In the text 
Fig. 3. Posterior distribution of the model parameters describing the annual modulation of the scintillation rate. The panels on the diagonal are the marginalised onedimensional posterior distributions. The dashed black lines are drawn at 0.16, 0.5, and 0.84 quartiles. The annotated text gives the bestfit parameters and ±1σ bounds implied by the posterior distributions. 

In the text 
Fig. 4. Scintillation rate determined through a Gaussian regression analysis (blue points with ±1σ errors); the maximum likelihood model (black curve) is overplotted. The annotated text gives the model parameters and the reduced χ^{2} of the fit. The phase at the vernal equinox is defined to be zero, with phase increasing with time. Magenta crosses show the inverse of the scintillation timescale measured directly from the autocorrelation function (ACF) of the light curves during epochs of sufficiently rapid variations. 

In the text 
Fig. 5. Left panel: line of sight and transverse distance of stars in the HIPPARCOS catalogue with respect to the sightline towards J1402+5347. The stars are colourcoded by their spectral type (legend in the annotated text). Right panel: same as the left panel, but the stars are plotted in the plane of the sky. The central black square shows J1402+5347. The solid and dashed blue lines mark the orientation of the major and minor axes of the scintels, respectively. 

In the text 
Fig. A.1. Light curves derived for April 9, 2019, of the variable source J1402+5347 and the steady comparison source J1403+5349, 9 arcmin west of J1402+5347, illustrating the accuracy of the method for constructing the light curves. The horizontal axis shows the time since the start of the observation. 

In the text 
Fig. B.1. Numerical ACF computation: mean subtracted light curve (top panel) and the ACF (bottom panel) for the April 9, 2019, light curve. 

In the text 
Fig. B.2. Bestfitting kernel: numerically computed ACF overplotted with the bestfit Gaussian kernel for the April 9, 2019, light curve. The decorrelation timescale of the Gaussian kernel is indicated with the dashed black line. The quasisinusoidal oscillations at long delay are not ideally captured by the analytic form of the kernel, However, the decay to the 1/evalue at short delays is well captured and is most relevant to estimating scintillation rates. 

In the text 
Fig. B.3. Gaussian process kernel parameters: posterior distribution of GP kernel parameters for the April 9, 2019, data set. The total variance K_{0} was fixed to equal the noisebias corrected variance of the light curve, and only the scale parameters c and d were estimated. The dashed lines are placed at the 16, 50, and 84^{th} quantiles of the marginalised onedimensional distributions. The fractional error in the estimation of parameter c is about 10%. 

In the text 
Fig. C.1. Same as Fig. 3, but for a twodimensional scintel model given in Eq. (C.1). 

In the text 
Fig. D.1. Separate light curves for the lower and upper halves of the observing band of the discovery observation on April 9, 2019. The blue line is the light curve for the frequency range 1311−1408 MHz, and the red line shows the range 1408−1505 MHz. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.