Issue 
A&A
Volume 628, August 2019



Article Number  A116  
Number of page(s)  8  
Section  Planets and planetary systems  
DOI  https://doi.org/10.1051/00046361/201731966  
Published online  14 August 2019 
Indications for transittiming variations in the exoNeptune HATP26b^{★}
^{1}
Stellar Astrophysics Centre, Department of Physics and Astronomy, Aarhus University,
Ny Munkegade 120,
8000 Aarhus C, Denmark
email: cessen@phys.au.dk
^{2}
Astronomical Observatory, Institute of Theoretical Physics and Astronomy, Vilnius University,
Sauletekio av. 3,
10257
Vilnius, Lithuania
^{3}
Institute of Theoretical Astrophysics, University of Oslo,
Postboks 1029 Blindern,
0315 Oslo,
Norway
^{4}
Facultad de Ciencias Astronómicas y Geofísicas, Universidad Nacional de La Plata,
Paseo del Bosque,
B1900FWA, La Plata, Argentina
^{5}
Instituto de Astrofísica de La Plata (CCTLa Plata, CONICETUNLP),
Paseo del Bosque, B1900FWA
La Plata, Argentina
^{6}
Complejo Astronómico El Leoncito (CONICETUNLPUNCUNSJ),
Av. España 1512 Sur, San Juan, Argentina
^{7}
Institut für Astrophysik, GeorgAugustUniversität Göttingen,
FriedrichHundPlatz 1,
37077
Göttingen, Germany
^{8}
Monash Centre for Astrophysics, School of Physics and Astronomy, Monash University,
Victoria 3800, Australia
^{9}
LeibnizInstitut für Astrophysik Potsdam,
An der Sternwarte 16,
14482
Potsdam, Germany
Received:
18
September
2017
Accepted:
15
April
2019
Upon its discovery, the lowdensity transiting Neptune HATP26b showed a 2.1σ detection drift in its spectroscopic data, while photometric data showed a weak curvature in the timing residuals, the confirmation of which required further followup observations. To investigate this suspected variability, we observed 11 primary transits of HATP26b between March, 2015, and July, 2018. For this, we used the 2.15 m Jorge Sahade Telescope placed in San Juan, Argentina, and the 1.2 m STELLA and the 2.5 m Nordic Optical Telescope, both located in the Canary Islands, Spain. To add to valuable information on the transmission spectrum of HATP26b, we focused our observations in the Rband only. To contrast the observed timing variability with possible stellar activity, we carried out a photometric followup of the host star over three years. We carried out a global fit to the data and determined the individual midtransit times focusing specifically on the light curves that showed complete transit coverage. Using bibliographic data corresponding to both ground and spacebased facilities, plus our new characterized midtransit times derived from partsperthousand precise photometry, we observed indications of transit timing variations in the system, with an amplitude of ~4 min and a periodicity of ~270 epochs. The photometric and spectroscopic followup observations of this system will be continued in order to rule out any aliasing effects caused by poor sampling and the longterm periodicity.
Key words: instrumentation: photometers / techniques: photometric / time / planets and satellites: fundamental parameters / stars: activity
The transit photometry (time, flux, error) and the long term monitoring in three bands are only available at the CDS via anonymous ftp to cdsarc.ustrasbg.fr (130.79.128.5) or via http://cdsarc.ustrasbg.fr/vizbin/qcat?J/A+A/628/A116
© ESO 2019
1 Introduction
From the first exoplanets ever detected around stars other than our Sun, the most successful exoplanet detection methods have been the radial velocity (RV; e.g., Butler et al. 2006) and the transit (Charbonneau et al. 2000, and onward) techniques. The main engine of this work, the transittiming variation (TTV) method, gained special relevance with the advent of the Kepler space telescope (Borucki et al. 2010; Koch et al. 2010). For this technique, once an exoplanet is detected via the transit method, the variations of the observed midtransit times permit subsequent detection and/or characterization of further transiting (Carter et al. 2012) and nontransiting (Barros et al. 2014) exoplanets. The technique is sensitive to planets with masses as low that of Earth (Agol et al. 2005), which would be extremely challenging to detect or characterize by means of other techniques. The TTV method requires sufficiently longbaseline and precise photometry, and good phase coverage. All these requirements are satisfied by Kepler data (see e.g., Mazeh et al. 2013, for a TTV characterization of hundreds of Kepler objects of interest, KOIs). As a consequence, surveys focused on TTVs from the ground have been mainly carried out photometrically followingup hot Jupiters with deep transits orbiting bright stars (e.g., Maciejewski et al. 2011; Fukui et al. 2011; von Essen et al. 2013). Nonetheless, none of them have revealed unquestionable detection of TTVs. On the contrary, TTVs in the Kepler era have revealed that multiple systems are not rare at all, and that most of the planets in these multiple systems are within the Super Earth/mini Neptune regime (see Holman et al. 2010; Lissauer et al. 2011; Cochran et al. 2011). Therefore, from Kepler results we can refocus our observing capabilities and boost our success rate from the ground by programming photometric followups of more suitable transiting systems, including the Neptune planets, instead of only hot Jupiters.
Between 2015 and 2018 our group carried out a photometric followup of three Neptunesized exoplanetary systems. The observations were mainly focused around primary transits. These systems presented either “hints” of multiplicity, or showed orbital and physical parameters similar to KOIs in multiple systems with largeamplitude TTVs. Here we present our efforts in the detection of TTVs in HATP26b, which is one of such aforementioned systems. HATP26b, a lowdensity Neptunemass planet transiting a Ktype star, and was discovered in 2011 by the HATNet consortia (Hartman et al. 2011). Once the spectral reconnaissance was carried out, additional highresolution spectra were acquired to better characterize the RV variations and the stellar parameters. The spectroscopic analysis revealed a chromospherically quiet star, along with a Neptunemass planet orbiting its host star every ~4.23 days. A combined analysis between photometric and spectroscopic data provided better constraints on the planetary mass and radius. Further analysis of spectroscopic data revealed a 2.1σ detection drift in its RVs. Nonetheless, the data were not sufficient to characterize the multiplicity of the system. Furthermore, Stevenson et al. (2016) observed a weak curvature in the timing residuals of HATP26b, but yet again without proper confirmation due to insufficient data. Although the ΔF = 0.53% transit depth makes observation of the transit challenging using groundbased facilities, the observed RV drift and curvature in the timing residuals makes HATP26b an interesting system for TTV studies. To confront the TTV detection with other possible sources of variability, we also carried out a photometric followup of the host star over three years.
Besides TTV analysis, several efforts were made to characterize the chemical constituents of the atmosphere of HATP26b. While Stevenson et al. (2016) found tentative evidence for water and a lack of potassium in transmission, Wakeford et al. (2017) measured the atmospheric heavy element content of the planet and characterized its atmosphere as primordial. The observations presented here focused solely on the Rband. Thus, besides the TTV discovery and characterization, in this work we also provide an accurate wavelengthdependent planettostar radii ratio around the Johnson–Cousins R band, 635 ± 100 nm, in order to contribute to its exoatmospheric characterization.
In this work, Sect. 2 details our observations and specifies the datareduction techniques. Section 3 shows the steps involved in the transit analysis and the derived model parameters. Section 4 shows our results on TTVs, and we finish with a brief discussion and our conclusions in Sect. 5.
2 Observations and data reduction
2.1 Observing sites and specifications of the collected data for transit photometry
Between March, 2015, and July, 2018, we observed HATP26 (R = 11.56, V = 11.76, Høg et al. 2000) during ten primary transits. Our observations were performed using the 2.15 m Jorge Sahade Telescope located at the Complejo Astronómico El Leoncito^{1} (CASLEO) in San Juan, Argentina, the 2.5 m Nordic Optical Telescope^{2} (NOT) located in La Palma, Spain, and the 1.2 m STELLA, located in Tenerife, Spain. The most relevant parameters of our observationsare summarized in Table 1. All the observations were carried out using a Johnson–Cousins R filter, to minimize the impact of our Earth’s atmosphere on the overall signaltonoise ratio (S/N) of our measurements, circumventing telluric linesand the strong absorption of stellar light caused by Rayleigh scattering. Also, several transits allowed us to accurately characterize the size of HATP26b within that wavelength range. To increase the photometric precision of our data, all telescopes were slightly defocused (Kjeldsen & Frandsen 1992; Southworth et al. 2009). In consequence, values of seeing taken from scienceframes are not representative of the quality of the sky at both sites. Table 1 shows then only intra night variability of the seeing. Exposure times were typically set between 10 and 60 s, and the photometricprecision ranged between 1.1 and 3.7 partsperthousand (ppt), in all cases below the transit depth (~5.3 ppt). Seven light curves show complete transit coverage. The rest show partial transit coverage mainly due to poor weather conditions.
During each observing night we obtained a set of 15 bias and between 10 and 15 twilight sky flats. Due to low exposure times and optimum cooling of chargecoupled devices, we did not take dark frames. In the case of CASLEO data we observed using a focal reducer. This produces an unvignetted field of view of 9 arcmin of diameter, allowing the simultaneous observation of HATP26 along with two reference stars of similar brightness, namely TYC 3204261 (V = 11.08, Høg et al. 2000) and TYC 32013761 (V = 12.5, Høg et al. 2000), also observed by STELLA. On the contrary, the size of the field of view of the NOT is about 7 arcmin^{2}. Therefore, we only observed TYC 3204261 and fainter reference stars simultaneously to HATP26.
2.2 Transit data reduction and preparation
For details on the data reduction, we refer the reader to the description of DIP^{2} OL (von Essen et al. 2018). In brief, all the science frames were bias and flatfield calibrated using the IRAF task ccdproc. Owing to the large availability of calibration frames, we corrected the science frames of a given observing night with the calibrations taken during that particular night only. To produce photometric light curves we used the IRAF task apphot. We measured fluxes inside apertures centered on the target and all available reference stars within the field of view of the telescopes. Their sizes were set as a fraction of the full width at half maximum (FWHM) computed and averaged per night. In particular, the apertures were set between 0.7 and 5 times the FWHM, divided into ten equally spaced steps. For each one of the apertures we chose three different sky rings, their widths being one, two, and three times the FWHM. The inner radii of the sky ring was fixed to five times the FWHM.
Next, we produced differential light curves for the target and each combination of reference stars by dividing the flux of the target by the average of a given combination of fluxes from the reference stars. For each one of the light curves we computed the standard deviation of the residuals, which were obtained by dividing the differential light curve by a highdegree, timedependent polynomial that was fitted to the whole data, including intransit points, minimizing the sum of the squares of the residuals. In particular, the degree of the polynomial was chosen to account for the shallow transit depth and most of the systematic noise. After visually inspecting all the light curves, it was chosen to be a septic degree and was fixed during all nights. The final combination of reference stars and aperture was chosen by minimizing the standard deviation of the residual light curves. This process was performed for each transit individually. With the time set in Julian dates and the differential fluxes defined, we finally changed the magnitude of the error bars provided by the apphot task so that their average magnitude was coincident with the standard deviation of the data. Besides time, flux, and errors, we also computed (x,y) centroid positions, integrated flat counts in the selected aperture around those (x,y) positions, the seeing, the airmass, and the background counts inside the selected sky ring, all of them as a function of time and for all the stars involved in the selected differential light curve.
Relevant parameters of our observations.
2.3 Datareduction of the longterm monitoring of HATP26
Stellar activity, and particularly stellar spots rotating with the star, can mimic TTVs (see Sect. 4.2). As part of the VAriability MOnitoring of exoplanet host Stars (VAMOS) project, we carried out a photometric followup of the host star, HATP26, over three years. To this end we observed with the widefield imager WIFSIP mounted at the 1.2 m twintelescope STELLA located in the Canary Islands, Spain (Strassmeier et al. 2004; Weber et al. 2012). Table 2 shows the main characteristics of the observations. The data were reduced using standard routines of ESOMIDAS. On average, around five stars were used to construct the differential light curve. For details on the data reduction steps, we refer the reader to Mallonn et al. (2015) and Mallonn & Strassmeier (2016).
3 Model parameters and transit analysis
3.1 Choice of detrending model and correlated noise treatment
To clean thedata from systematic errors we constructed an initial detrending model, DM, with its terms represented by a linear combination of seeing, S, airmass, χ, centroid positions of all the stars involved in the creation of the differential light curves (x,y), and their respective integrated flat counts and background counts around these positions, henceforth, the detrending components: (1)
Here, N is the total number of reference stars and N + 1 includes the components of the target star. X_{i} and Y_{i} correspond to the (x,y) centroid positions, FC_{i} to the integrated flat counts, and BG_{i} to the background counts. The coefficients of the detrending model are a_{0}, a_{1}, a_{2}, and d_{i}, e_{i} and f_{i}, with i = 1, N + 1. Besides this linear combination, we also considered timedependent polynomials with degrees ranging from one to three. The detrending strategy is fully described in von Essen et al. (2018). In an iterative process we defined a joint model described as the transit model (see Sect. 3.2 for details) multiplied by the detrending model. During each iteration, where we evaluated the different timedependent polynomials and the different combinations of detrending components, we computed reduced χ^{2}, the Bayesian information criterion (BIC; see e.g., von Essen et al. 2017, for details on its use) and the weighted standard deviation between the transit photometry and the combined model. In all these statistical indicators, the number of fitting parameters is taken into consideration. Thus, we made sure that the data were not being overfitted by an unnecessarily large detrending model. After averaging the three statistical indicators, we chose as final detrending model (and thus, the final number of detrending components or the final degree for the timedependent polynomial) the one that minimized this averaged statistic. Usually, the airmass, the (x,y) centroid positions, and the integrated flat counts of all the stars involved in the construction of the differential light curves were part of the chosen detrending model. For the NOT data, which counts with an excellent tracking system, a secondorder, timedependent polynomial was usually sufficient to account for instrumental systematics. Once the detrending model was defined, we treated correlated noise as specified in von Essen et al. (2013). Briefly, following Carter & Winn (2009) we produced residual light curves, subtracting the primary transit model from the data (Mandel & Agol 2002). As orbital parameters we used bibliographic values, listed in this work in the first column of Table 3. We then divided each light curve into M bins of equal duration, each bin holding an averaged number of N data points.An averaged N accounts for unequally spaced data points. If the data are free of correlated noise, then the noise within the residual light curves should follow the expectation of independent random numbers , (2)
Here, σ_{1} is the sample variance of the unbinned data and σ_{N} is the sample variance (or RMS) of the binned data: (3)
where the mean value of the residuals per bin is given by , and is the mean value of the means. When correlated noise cannot be neglected, each σ_{N} differs by some factor β_{N} from their expectation . We computed β averaging the values of β_{N} obtained over time bins close to the duration of the transit ingress (or, equivalently, egress) of HATP26b, which is estimated to be ~15 min. To estimate β, we used time bins of 0.8, 0.9, 1, 1.1, and 1.2 times the duration of ingress. Finally, we enlarged the magnitude of the error bars by β (Pont et al. 2006). On average, CASLEO data presented slightly more correlated noise (β ~ 2.1) compared to the NOT (β ~ 1.5). This is explained by guiding the stars (see column eight of the Table 1 for a measurement of the amplitude of variability of the centroid position of the stars). While the NOT has a stable guiding system, CASLEO lacks guiding and the telescope slightly tracks in the eastwest direction, affecting the quality of the photometry.
Main characteristics of the observations performed with STELLA.
Bibliographic and derived transit parameters of HATP26b.
3.2 Transit fitting
With the ten light curves fully constructed, we converted the time stamps from Julian dates to Barycentric Julian dates, BJD_{TDB}, using the web tool provided by Eastman et al. (2010). The first step was to determine the expectation values of the model parameters of HATP26b. To this end, we carried out a Markovchain Monte Carlo (MCMC) global fit using the PyAstronomy^{3} packages (Patil et al. 2010; Jones et al. 2001). As transit model we used the one provided by Mandel & Agol (2002), along with a quadratic limbdarkening law. In particular, for the Mandel & Agol (2002) transit model the fitting parameters are the inclination in degrees, i, the semimajor axis in stellar radius, a∕R_{S}, the orbital period in days, Per, the planettostar radii ratio, R_{P} /R_{S}, the midtransit time in BJD_{TDB}, T_{0}, and the linear and quadratic limbdarkening coefficients, u_{1} and u_{2}.
The quadratic limbdarkening coefficients were taken from Claret (2000) for the Johnson–Cousins R filter and stellar parameters closely matching those of HATP26 (T_{eff} = 5079 ± 88 K, [Fe/H] = −0.04 ± 0.08, and log (g) = 4.56 ± 0.06, Hartman et al. 2011). Thus, the values for the limbdarkening coefficients used in this work correspond to the stellar parameters T_{eff} = 5000 K, [Fe/H] = 0 and log (g) = 4.5, and are listed in Table 3. In this work, quadratic limbdarkening coefficients are considered as fixed values. Nonetheless, to test if the adopted values for the limbdarkening coefficients have an impact on the measured midtransit times, we considered the linear and quadratic limbdarkening values as fitting parameters as well. After carrying out a global fit including them in the fitting process, our posterior distributions for the parameters revealed values inconsistent with the spectral classification of the host star. Although our data present a photometric precision at the partsperthousand level, the transit light curves are not precise enough to realistically fit for the limbdarkening coefficients. Nonetheless, midtransit times computed fitting for the limbdarkening coefficients do not show a significant offset with respect to those computed to fix them, but they do show slightly larger error bars.
To fit the data with our combined model, we began our analysis considering the transit parameters listed in Hartman et al. (2011), Stevenson et al. (2016), and Wakeford et al. (2017) as starting values. The adopted values are shown in the second column of Table 3. We chose uniform priors for the orbital period, and Gaussian priors for the semimajor axis, the inclination, and the planettostar radii ratio. To be conservative and to avoid misleading the convergence of the fitting algorithm, as standard deviation for the Gaussian priors we chose three times the magnitude of the errors determined by the latter mentioned authors. The midtransit time determined by Stevenson et al. (2016) has already five decimals of precision, which translates into an accuracy of ~2 s. Since this precision is significantly below our data cadence, throughout this work the reference midtransit time is consider as fixed.
Once the priors were set, we iterated the fitting process 1.5 ×10^{5} times. At each MCMC step not only are the transit parameters changed, but also the detrending coefficients are computed with the previously mentioned inversion techniques. Both transit parameters and detrending coefficients are saved at each MCMC step. When the total number of iterations was reached, we burned the initial 50 000 samples. From the mean values and standard deviations of the samples we computed the expectation values of the model parameters and their corresponding errors. The same process was followed to determine the bestfit detrending coefficients. Throughout this work we provide the errors of the transit parameters at 1σ level. As a consistency check we analyzed the posterior distributions to ensure convergence of the chains. To this end, we divided the chains in three equally large fractions, we computed the model parameters for each one of the subchains, and we checked that they were consistent within 1σ errors. Our bestfit values are listed in the third column of Table 3. To avoid visual contamination, not all the values for the transit parameters derived by other authors are displayed in Table 3. Nonetheless, the transit parameters presented in this work are in good agreement with all the bibliographic values. In particular, the value for the planettostar radii ratio follows the trend observed by Wakeford et al. (2017), and not the value reported by Hartman et al. (2011). Figure 1 shows one transit light curve and both model components in isolation.
After producing a global fit to the complete data set (see Sect. 3.2), to characterize possible TTVs in the system we proceeded to fit the transits individually. Figure 2 shows the 11 transits, 7 of them corresponding to the JS telescope(red circles), 1 to STELLA (blue circles) and the remaining 4 to the NOT (green squares). To provide TTVs that are as realistic as possible, for this analysis we fitted each individual transit. However, for a better assessment of the TTVs, complete and incomplete transits are clearly differentiated.
To compute the midtransit times of the individual light curves, rather than considering our bestfit model parameters as fixed, we used our global bestfit values and errors for the mean and standard deviation, respectively, for the Gaussian priors. For each individual midtransit time we chose uniform priors. Since the determination of the individual midtransit times requires the analysis of individual light curves, the orbital period was considered as fixed to the value derived from the global fit. This will properly propagate the uncertainties of all the transit parameters into the determination of the midtransit times and thus will provide more realistic uncertainties. Equivalently to the global fit, we iterated and burned the same number of chains per light curve, and we computed the individual midtransit times in the same fashion as specified in Sect. 3.2. For each light curve we also computed the transit parameters and checked that they were consistent with the global ones, always at 1σ level. The first section of Table 4 shows the individual midtransit times from the literature and their respective uncertainties, while the second section shows the midtransit times and the uncertainties derived in this work.
Fig. 1
Detrending process. The figure shows the evolution of the transit photometry during the fitting procedure. Raw data are plotted in red circles. The detrending model is plotted as a black continuous line. Artificially shifted, the detrended data are plotted as blue squares, along with primary transit models whose transit parameters are determined by the posterior distributions of the MCMC chains. 
4 Analysis and results
4.1 Transit timing variations
We computed the TTVs for the HATP26b system by subtracting observed midtransit times from midtransit times obtained from a constant (bestfit) orbital period. The derived values can be found in Fig. 3and the lower part of Table 4. We note that in the figure the zeroth epoch does not lie over the abscissas axis. This is due to the improvement of its value between Hartman et al. (2011) and Wakeford et al. (2017).
In the absence of any timing variations we expect no significant deviations of the derived O–C values from zero. We used a test to test the null hypothesis of O–C = 0. For 20 degrees of freedom the value raises up to = 177, with a pvalue of P < 1 × 10^{−5}. For this, the significance level is chosen to be 1%. This simple statistical analysis, , and its corresponding pvalue indicate that the midtransit times of HATP26b are inconsistent with a constant period. Visual inspection of the O–C diagram for HATP26b suggests that the TTVs are strongly dependent on the last three data points. To test these assumptions we excluded the three last epochs and recomputed our statistics: = 102, P < 1 × 10^{−5}.
For completeness, assuming TTVs with a sinusoidal variability we fitted a simple sinusoidal function with the following expression to our timing residuals: (4)
where E corresponds to the epoch, A_{TTVs} corresponds to the semiamplitude of the timing residuals, Per_{TTVs} to the periodicity, and ϕ_{TTVs} to the phase. From our analysis,we obtain A_{TTVs} ~ 2.1 min, while the bestfit periodicity would correspond to around 270 epochs. This is equivalent to almost 1100 days. Comparing the semiamplitude with the average error bars for the timings we find A_{TTVs} ~ 3 × ɛ, ɛ being the averaged value of the timing errors. Our bestfit sinusoidal function is plotted in Fig. 3 as a black dashed line. Mainly due to possible aliasing effects, and due to the complications that singletransiting planets with TTVs convey, we believe any analysis on the TTVs reported in this work is subject to strong speculation. More spectroscopic and photometric data are required to characterize the body causing these TTVs.
4.2 Transit timing variations induced by stellar activity
Due to the high photometric quality provided by spacebased observations such as CoRoT (Auvergne et al. 2009) and Kepler, stellar magnetic activity and its impact on transit light curves have been studied in great detail. Dark spots and bright feculae on the stellar photosphere move as the star rotates, producing a timedependent variation of the stellar flux. The currently obtainable photometric precision is so high that the small imprints of occulted spots on transit data have been used to characterize stellar surfacebrightness profiles and spot migration and evolution (see e.g., Carter et al. 2011; SanchisOjeda et al. 2011; SanchisOjeda & Winn 2011; Huber et al. 2010). When transit fitting is performed, an incorrect treatment of these “bumps” can have a significant impact on the computation of planetary and stellar parameters (see e.g., Czesla et al. 2009; Lanza 2011).
Occulted and unocculted spots can asymmetrically modify the shape of the transit light curves, and thus affect the true locations of the midtransit times. Actually, the deformations that stellar activity produces on primary transits have been studied in detail and pinpointed in some cases as a misleading identification of TTVs (see e.g., Rabus et al. 2009). For example, Maciejewski et al. (2011) found TTVs in the WASP10 system with an amplitude of a few minutes. These latter authors attributed the variability to the gravitational interaction between the transiting planet and a second body of a tenth of a Jupiter mass close to a 5:3 meanmotion resonance. Although Barros et al. (2013) did not find a linear ephemeris as a statistically good fit to the midtransit times of WASP10b, they showed that the observed variability could be produced by, for example, spot occultation features.
To quantify stellar activity in our system and support (or overrule) our TTV detection, we carried out a photometric followup of HATP26 over three years (see e.g., Mallonn et al. 2015; Mallonn & Strassmeier 2016, Fig. 4 and our Table 2). The data, taken in the Johnson–Cousins B, V, and I filters, show a maximum scatter of 2.3 ppt. To search for any significant periodicity contained in the data we computed Lomb–Scargle periodograms (Lomb 1976; Scargle 1982; Zechmeister & Kürster 2009), finding featureless periodograms within each observing season and within the whole observing run. As a consequence, the photometric followup of the host star seems to show no evidence of spot modulation. Nonetheless, the characterization presented in this work is as good as the data are. In other words, if spotinduced modulation were to exist, it would be well contained within the ~2 ppm limit.
Hartman et al. (2011) characterized the chromospheric activity of the star derived from the flux in the cores of the Ca II H and K lines, the socalled S index, S_{HK} = 0.182 ± 0.004. Comparing this value to the relation between the S index and the stellar brightness variation found by Karoff et al. (2016) (their Fig. 5), observational evidence should place the photometric variability of HATP26 at around 1 ppt. Assuming that spot modulation around and below this limit actually exists, and that it might have an impact on the TTVs, Ioannidis et al. (2016) found that the maximum amplitude of TTVs generated by spot crossing events strongly depends on transit duration (among other things). In the case of HATP26b (~143 min), TTVs caused by spots should have a maximum amplitude of ~1 min. This is significantly below the TTV amplitude detected in this work. Therefore, the data presented here seem to support a scenario where TTVs are caused by gravitational interactions instead of being modulated by the potential activity of the star.
Fig. 2
11 transit light curves collected in this work. The transits observed with the Argentinian 2.15 m telescope, CASLEO, are plotted as red circles. The transits taken with the Scandinavian 2.5 m NOT are shown as green squares. The transit observed with STELLA is shown in blue circles. Calender dates are displayed, and asterisks indicate the transits that are complete. Vertical dashed lines indicate the four contact times to guide the reader. Transits have been folded using our bestfit orbital period. In consequence, TTVs can be appreciated by comparing ingress or egress to the vertical dashed lines. 
Fitted instant of minima with 1σ errors and the O–C data points.
Fig. 3
Observed minus calculated (O–C) midtransit times for HATP26b in minutes. Red circular points show the values obtained from Hartman et al. (2011), Stevenson et al. (2016), and Wakeford et al. (2017). Black filled squares correspond to the values derived in this work considering complete transits, while empty squares are those derived from incomplete transits. The O–C diagram was computed using the bestfit orbital period derived in this work, specified in Table 3. The dashed black line corresponds to our bestfit sinusoidal variability. 
Fig. 4
Photometric followup of HATP26 using the STELLA photometer. Red circles correspond to the I band, green squares to the V band, and blue triangles to the B band. The fluxes have been artificially shifted, and the time axis has been compressed for improved visualization. Dashed black lines show the averaged standard deviation. 
4.3 Transit timing variations caused by systematic effects
Besides stellar activity, TTVs can be caused by systematic effects not properly accounted for. This has a special relevance when lowamplitude TTVs are being reviewed. A manifestation of this effect is the derivation of underestimated error bars for the midtransit times. While in this work we have taken special care over the computation of uncertainties, we can only trust that the values reported by other authors have received similar treatment. Nonetheless, as a consistency check we enlarged the error bars of all the observed midtransit times by a factor and recomputed for each one of the artificially enlarged O–C diagrams. A factor of 13 was required for to be consistent with 1, equivalently, with no TTVs. This marginally relatively large value seems to support our TTV detection.
5 Discussion and conclusions
Since its discovery, indications of variability in both spectroscopic and photometric data have flagged HATP26b as an interesting target for additional photometric followup. For the last three years our group has been collecting primary transit data of partsperthousand photometric precision, necessary to detect the shallow transits that the planet produces while crossing its host star every ~4.23 days. In this study we reduced all the new data in a homogeneous way, updated and improved the transit parameters and ephemeris, and detected significant variability in the timing residuals. Furthermore, we followed the host star over several years to characterize its activity and, if observed, its impact on midtransit times. However, light curves taken in three different bands over three years revealed no spot modulation within the precision limit of the data, devaluing the scenario of spotinduced TTVs. Due to the complexity of the analysis of singletransiting systems presenting TTVs, it is hard to make any robust conclusions on the characteristics of the perturbing body. Nonetheless, we understand these results to be relevant and to motivate future photometric and spectroscopic followup campaigns that will contribute to disclosing the origin of the observed variability.
Acknowledgements
Visiting Astronomer, Complejo Astronómico El Leoncito operated under agreement between the Consejo Nacional de Investigaciones Científicas y Técnicas de la República Argentina and the National Universities of La Plata, Córdoba and San Juan. C.v.E. acknowledges funding for the Stellar Astrophysics Centre, provided by The Danish National Research Foundation (Grant DNRF106), and support from the European Social Fund via the Lithuanian Science Council grant No. 09.3.3LMTK712010103. This study is based on observations made with the NOT, operated by the Nordic Optical Telescope Scientific Association at the Observatorio del Roque de los Muchachos, La Palma, Spain, of the Instituto de Astrofisica de Canarias. We acknowledge support from the Research Council of Norway’s grant 188910 to finance service observing at the NOT. S.W. acknowledges support for International Team 265 (“Magnetic Activity of Mtype Dwarf Stars and the Influence on Habitable Extrasolar Planets”) funded by the International Space Science Institute (ISSI) in Bern, Switzerland, and support by the Research Council of Norway through its Centres of Excellence scheme, project number 262622. V.P. wishes to acknowledge the funding that was provided in part by a J. William Fulbright Grant and the sponsorship of GeorgAugustUniversitat in Gottingen, Germany. We thank the referee for comments which contributed to improving this work. This work made use of PyAstronomy^{4}.
References
 Agol, E., Steffen, J., Sari, R., & Clarkson, W. 2005, MNRAS, 359, 567 [NASA ADS] [CrossRef] [Google Scholar]
 Auvergne, M., Bodin, P., Boisnard, L., et al. 2009, A&A, 506, 411 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Barros, S. C. C., Boué, G., Gibson, N. P., et al. 2013, MNRAS, 430, 3032 [NASA ADS] [CrossRef] [Google Scholar]
 Barros, S. C. C., Díaz, R. F., Santerne, A., et al. 2014, A&A, 561, L1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Borucki, W. J., Koch, D., Basri, G., et al. 2010, Science, 327, 977 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Butler, R. P., Wright, J. T., Marcy, G. W., et al. 2006, ApJ, 646, 505 [NASA ADS] [CrossRef] [Google Scholar]
 Carter, J. A., & Winn, J. N. 2009, ApJ, 704, 51 [NASA ADS] [CrossRef] [Google Scholar]
 Carter, J. A., Winn, J. N., Holman, M. J., et al. 2011, ApJ, 730, 82 [NASA ADS] [CrossRef] [Google Scholar]
 Carter, J. A., Agol, E., Chaplin, W. J., et al. 2012, Science, 337, 556 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Charbonneau, D., Brown, T. M., Latham, D. W., & Mayor, M. 2000, ApJ, 529, L45 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Claret, A. 2000, A&A, 363, 1081 [NASA ADS] [Google Scholar]
 Cochran, W. D., Fabrycky, D. C., Torres, G., et al. 2011, ApJS, 197, 7 [NASA ADS] [CrossRef] [Google Scholar]
 Czesla, S., Huber, K. F., Wolter, U., Schröter, S., & Schmitt, J. H. M. M. 2009, A&A, 505, 1277 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Eastman, J., Siverd, R., & Gaudi, B. S. 2010, PASP, 122, 935 [NASA ADS] [CrossRef] [Google Scholar]
 Fukui, A., Narita, N., Tristram, P. J., et al. 2011, PASJ, 63, 287 [NASA ADS] [Google Scholar]
 Hartman, J. D., Bakos, G. Á., Kipping, D. M., et al. 2011, ApJ, 728, 138 [NASA ADS] [CrossRef] [Google Scholar]
 Høg, E., Fabricius, C., Makarov, V. V., et al. 2000, A&A, 355, L27 [Google Scholar]
 Holman, M. J., Fabrycky, D. C., Ragozzine, D., et al. 2010, Science, 330, 51 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Huber, K. F., Czesla, S., Wolter, U., & Schmitt, J. H. M. M. 2010, A&A, 514, A39 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Ioannidis, P., Huber, K. F., & Schmitt, J. H. M. M. 2016, A&A, 585, A72 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Jones, E., Oliphant, T., Peterson, P., et al. 2001, SciPy: Open source scientific tools for Python, http://www.scipy.org [Google Scholar]
 Karoff, C., Knudsen, M. F., De Cat, P., et al. 2016, Nat. Commun., 7, 11058 [NASA ADS] [CrossRef] [Google Scholar]
 Kjeldsen, H., & Frandsen, S. 1992, PASP, 104, 413 [NASA ADS] [CrossRef] [Google Scholar]
 Koch, D. G., Borucki, W. J., Basri, G., et al. 2010, ApJ, 713, L79 [NASA ADS] [CrossRef] [Google Scholar]
 Lanza, A. F. 2011, in Physics of Sun and Star Spots, eds. D. Prasad Choudhary & K. G. Strassmeier, IAU Symp., 273, 89 [NASA ADS] [Google Scholar]
 Lissauer, J. J., Fabrycky, D. C., Ford, E. B., et al. 2011, Nature, 470, 53 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Lomb, N. R. 1976, Ap&SS, 39, 447 [NASA ADS] [CrossRef] [Google Scholar]
 Maciejewski, G., Dimitrov, D., Neuhäuser, R., et al. 2011, MNRAS, 411, 1204 [NASA ADS] [CrossRef] [Google Scholar]
 Mallonn, M., & Strassmeier, K. G. 2016, A&A, 590, A100 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Mallonn, M., von Essen, C., Weingrill, J., et al. 2015, A&A, 580, A60 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Mandel, K., & Agol, E. 2002, ApJ, 580, L171 [NASA ADS] [CrossRef] [Google Scholar]
 Mazeh, T., Nachmani, G., Holczer, T., et al. 2013, ApJS, 208, 16 [NASA ADS] [CrossRef] [Google Scholar]
 Patil, A., Huard, D., & Fonnesbeck, C. J. 2010, J. Stat. Softw., 35, 1 [CrossRef] [EDP Sciences] [Google Scholar]
 Pont, F., Zucker, S., & Queloz, D. 2006, MNRAS, 373, 231 [NASA ADS] [CrossRef] [Google Scholar]
 Rabus, M., Alonso, R., Belmonte, J. A., et al. 2009, A&A, 494, 391 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 SanchisOjeda, R., & Winn, J. N. 2011, ApJ, 743, 61 [NASA ADS] [CrossRef] [Google Scholar]
 SanchisOjeda, R., Winn, J. N., Holman, M. J., et al. 2011, ApJ, 733, 127 [NASA ADS] [CrossRef] [Google Scholar]
 Scargle, J. D. 1982, ApJ, 263, 835 [NASA ADS] [CrossRef] [Google Scholar]
 Southworth, J., Hinse, T. C., Jørgensen, U. G., et al. 2009, MNRAS, 396, 1023 [NASA ADS] [CrossRef] [Google Scholar]
 Stevenson, K. B., Bean, J. L., Seifahrt, A., et al. 2016, ApJ, 817, 141 [NASA ADS] [CrossRef] [Google Scholar]
 Strassmeier, K. G., Granzer, T., Weber, M., et al. 2004, Astron. Nachr., 325, 527 [NASA ADS] [CrossRef] [Google Scholar]
 von Essen, C., Schröter, S., Agol, E., & Schmitt, J. H. M. M. 2013, A&A, 555, A92 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 von Essen, C., Cellone, S., Mallonn, M., et al. 2017, A&A, 603, A20 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 von Essen, C., Ofir, A., Dreizler, S., et al. 2018, A&A, 615, A79 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Wakeford, H. R., Sing, D. K., Kataria, T., et al. 2017, Science, 356, 628 [NASA ADS] [CrossRef] [Google Scholar]
 Weber, M., Granzer, T., & Strassmeier, K. G. 2012, in Software and Cyberinfrastructure for Astronomy II, Proc. SPIE, 8451, 84510K [CrossRef] [Google Scholar]
 Zechmeister, M., & Kürster, M. 2009, A&A, 496, 577 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Visiting Astronomer, Complejo Astronómico El Leoncito operated under agreement between the Consejo Nacional de Investigaciones Científicas y Técnicas de la República Argentina and the National Universities of La Plata, Córdoba and San Juan; www.casleo.gov.ar
All Tables
All Figures
Fig. 1
Detrending process. The figure shows the evolution of the transit photometry during the fitting procedure. Raw data are plotted in red circles. The detrending model is plotted as a black continuous line. Artificially shifted, the detrended data are plotted as blue squares, along with primary transit models whose transit parameters are determined by the posterior distributions of the MCMC chains. 

In the text 
Fig. 2
11 transit light curves collected in this work. The transits observed with the Argentinian 2.15 m telescope, CASLEO, are plotted as red circles. The transits taken with the Scandinavian 2.5 m NOT are shown as green squares. The transit observed with STELLA is shown in blue circles. Calender dates are displayed, and asterisks indicate the transits that are complete. Vertical dashed lines indicate the four contact times to guide the reader. Transits have been folded using our bestfit orbital period. In consequence, TTVs can be appreciated by comparing ingress or egress to the vertical dashed lines. 

In the text 
Fig. 3
Observed minus calculated (O–C) midtransit times for HATP26b in minutes. Red circular points show the values obtained from Hartman et al. (2011), Stevenson et al. (2016), and Wakeford et al. (2017). Black filled squares correspond to the values derived in this work considering complete transits, while empty squares are those derived from incomplete transits. The O–C diagram was computed using the bestfit orbital period derived in this work, specified in Table 3. The dashed black line corresponds to our bestfit sinusoidal variability. 

In the text 
Fig. 4
Photometric followup of HATP26 using the STELLA photometer. Red circles correspond to the I band, green squares to the V band, and blue triangles to the B band. The fluxes have been artificially shifted, and the time axis has been compressed for improved visualization. Dashed black lines show the averaged standard deviation. 

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.