LSPM J1112+7626: DETECTION OF A 41 DAY M-DWARF ECLIPSING BINARY FROM THE MEARTH TRANSIT SURVEY

LSPM J1112+7626 was targeted as part of routine operation of the MEarth transit survey (Nutzman & Charbonneau 2008; Irwin et al. 2009a) during the 2009–2010 season, with the first observations taken on UT 2009 December 2. Exposure times for the LSPM J1112+7626 field were 27 s, taking three exposures per visit (the total exposure times are tailored for each target to achieve sensitivity to a particular planet size for the assumed stellar radius), with visits at a cadence of 20 minutes.

Eclipses were first detected in LSPM J1112+7626 on UT 2010 April 3 (this was a primary eclipse), and subsequently confirmed on UT 2010 April 28 by observing the majority of a secondary eclipse egress. Both eclipses were detected by an experimental automated real-time eclipse detection and follow-up system we began to trial during the 2009–2010 season (this will be described in detail in a future publication when it has matured). Confident that the object was an EB, but still with an ambiguous orbital period, we immediately commenced radial velocity monitoring. The radial velocities ultimately determined which of the possible orbital periods was correct, in conjunction with a further eclipse (which was also a secondary) observed on UT 2010 June 8 with MEarth and the Clay telescope (see Section 2.2).

By the time the orbital ephemeris was well known, the observing season on LSPM J1112+7626 from Arizona was ending, and the weather had begun to deteriorate due to the arrival of the monsoon. We therefore took advantage of the high declination of our target and commenced an observing program to target the secondary eclipses using the Harvard University 0.4 m Landon T. Clay telescope, which is located on the roof of the Science Center on the Harvard campus near Harvard Square, Cambridge, Massachusetts. This instrument is normally used for undergraduate teaching. It is equipped with an e2v CCD47-10 1024 × 1024 pixel back-illuminated CCD packaged by Apogee Instruments (model AP47), which was operated binned 2 × 2 for a final pixel scale of ≈159 pixel⁻¹ to reduce overheads.

The location of the telescope is extremely light polluted, and combined with the red colors of our target, we decided to use the I-band filter to mitigate the effects of sky background and atmospheric extinction. The filter available was built to the prescription of Bessell (1990b) using Schott RG9 glass to define the bandpass, and shall hereafter be referred to as I_Bessell. While this filter was intended to approximate the Cousins I passband (hereafter I_C), the combination of the glass filter and CCD quantum efficiency yields a system response with significant sensitivity at very red wavelengths not present in the original Cousins passband, and thus is a poor approximation to the Cousins filter for very red stars such as M-dwarfs. We shall return to this issue in Section 4.4.

Data were successfully gathered using the Clay telescope for the secondary eclipses on UT 2010 June 8, October 9, November 19, and December 30. The entire eclipse was visible for all but the June 8 event, although varying fractions of these were lost due to clouds. Exposure times were adjusted to keep the peak counts in the target below 30,000 to avoid saturation and ranged from 30 to 60 s. The typical cadence was 50 s and the FWHM of the stellar images was 2–3 binned pixels. No guiding was available, so the target was re-centered manually approximately every 20 exposures to attempt to keep the drift below 5 pixels.

Data were reduced and light curves produced using standard procedures for CCD data, following the method outlined in Irwin et al. (2007). Two-dimensional bias frames were subtracted as there was some bias structure evident, and we used the overscan region to correct for any bias level variations. Dark frames of exposure matched to the target frames were subtracted to correct for the (significant) dark current and hot pixels at the −20°C operating temperature of the CCD, and the data were flat-fielded using twilight flat-fields in the usual way. No defringing or absolute calibration of the photometry were attempted. Fringing was clearly visible in the images, but at a low amplitude (3% peak-to-peak) and should not have a significant effect on the photometry due to the stabilized telescope pointing.

We present the full light curve data set used in this paper in Table 1, and times of minimum light for the eclipses where a clear minimum was seen in Table 2.

Using the timing of the single visit within the primary eclipse of UT 2010 April 3 from MEarth (see Section 2.1) combined with modeling of the secondary eclipses (Section 2.2) and radial velocities (Section 3) it was possible to predict the primary eclipse times with reasonable accuracy. During 2010–2011 full eclipses were observable at reasonable airmass (<2.0) from eastern Europe, Scandinavia, and western Asia, with the high declination of the target favoring northerly latitudes during winter to maximize dark hours at reasonable airmass.

The primary eclipse was observed using the 40 cm telescope at Hankasalmi Observatory, Finland on UT 2011 January 15, starting approximately ten minutes before first contact, and continuing until well after last contact. Data were taken using a similar I_Bessell filter to the secondary eclipse observations (see Section 2.2) with a front-illuminated Kodak KAF-1001E CCD packaged by Santa Barbara Instrument Group (model STL-1001E). This system response should be very similar to the one used for the secondary eclipses, but nonetheless, any bandpass mismatch is a potential concern for the accuracy of the photometric solution when combining data taken with different instruments, and will be discussed in more detail in Section 4.4.

Exposure times of 60 s were used, yielding a cadence of approximately two minutes, and the target was re-centered every ten exposures to keep the drift below ≈50 pixels. Data were reduced using two-dimensional bias frames, and master dark and flat-field frames. Light curves were then derived using the same methods and software as for the other data sets.

Following the discovery of atmospheric water vapor induced systematic effects in the MEarth photometry (see Irwin et al. 2011), a customized I-band filter was procured for the start of the 2010–2011 season, with a bandpass designed to eliminate the strong telluric water vapor absorption at wavelengths > 900 nm, by using an 895 nm (50% transmittance point) interference cutoff on the same RG715 glass substrate as the original filter (where the long-pass glass filter defines the blue end of the bandpass). This was predicted to reduce the effect to largely insignificant levels, and the resulting system response approximates the Cousins I bandpass (as tabulated by Bessell 1990b; see later in this section).

Observations of LSPM J1112+7626 using this new filter commenced on UT 2010 October 29 and continued until UT 2011 May 18. Exposures of 2 × 97 s were used at 20 minute cadence for times out of eclipse. Five thousand two hundred and eighty-three useful data points were obtained at airmasses <2, after discarding 400 bad frames (pointing errors, severe light loss due to clouds, and large scatter in the zero-point solution for the differential photometry).

Partial secondary eclipses were obtained with the same MEarth telescope on UT 2010 November 19, UT 2011 February 9, and UT 2011 May 2. Due to known systematic effects caused by persistence in the MEarth CCDs, leading to offsets in magnitude depending on the cadence (see Berta et al. 2011), and to remove any residual from the out-of-eclipse variations, we allow these high-cadence observations to have different photometric zero points than the out-of-eclipse data (and each other) in the final light curve analysis.

The out-of-eclipse data (excluding the high-cadence eclipse observations) were searched for periodic variations using the method outlined in Irwin et al. (2011). Correction for the "common-mode" effect was still performed, as we discovered during the season that replacing the filters had in fact increased the size (and changed the sign) of the systematic effect rather than eliminating it. This problem is still under investigation, but we suspect it is due to humidity and/or temperature dependence of the filter bandpass, an effect also seen in the original Sloan Digital Sky Survey (SDSS) Photometric Telescope filters (Doi et al. 2010).

After correcting for the "common-mode" effect, the MEarth data show a strong near-sinusoidal modulation with a 65 day period (see Figure 1; the data are shown in Section 4) and a peak-to-peak amplitude of approximately 0.02 mag. Marginal evidence for a similar periodicity is found in the 2009–2010 data, although the phase is not consistent and the amplitude appears to be much smaller. It is unclear if the amplitude change results from the differing filter systems, but the change in phase (if real) may be indicative of evolution of the photospheric spot patterns giving rise to the modulations.

Zoom In Zoom Out Reset image size

An improvement resulting from the new MEarth filter is a greatly reduced color term when standardizing the photometry onto the Cousins system. The scatter in this relation also appears to be significantly lower, which may also result in part from improvements made to our flat-fielding strategy during the 2010–2011 season. The MEarth observing software automatically obtains measurements of equatorial photometric standard star fields from the catalog of Landolt (1992) every night, and these are used to derive photometric zero points and calibrated photometry. The typical scatter in these zero-point solutions on a photometric night is 0.02 mag. We find a typical atmospheric extinction of 0.05 mag/airmass from MEarth data, which has been assumed henceforth.

A significant issue with the equatorial standard fields is they contain very few red stars suitable for defining the red end of the photometric response, and so color equations derived only from these equatorial fields may be significantly in error when attempting to standardize M-dwarf photometry. Thankfully, MEarth serendipitously targeted a number of M-dwarfs with photometry available on the Cousins system from Bessel (1990a). We selected all such objects in the 2010–2011 MEarth database with data available on photometric nights for calibration. Two objects with clearly unreliable MEarth data were discarded, leaving eight objects. The following color equation was derived from these stars, which span I_C − K_2MASS colors of 2.2–4.2:

$$\begin{eqnarray*} I_C&=& I_{\rm MEarth}- (0.085 \pm 0.019) \\ &&+ (0.0204 \pm 0.0064)\ (I_{\rm MEarth}-K_{\rm 2MASS}). \end{eqnarray*}$$

In the case of LSPM J1112+7626, we estimate I_C = 12.14 ± 0.05 using this relation, where we have inflated the uncertainty to account for variability and any errors introduced by the varying MEarth bandpass as discussed above. The observed I_C − J, I_C − H, and I_C − K colors are later used to infer effective temperatures and luminosities for the members of the binary, and thus the distance to the system.

In Table 3, we summarize the photometric and astrometric system properties gathered from literature data and our photometry. This star lies within the area covered by the Sloan Extension for Galactic Understanding and Evolution stripe 1260 photometry, so we retrieved the point-spread function (PSF) magnitudes from Data Release 7 (Abazajian et al. 2009). This object is flagged as saturated in i-band (the image also shows a very unusual PSF and may be corrupted), and the u-band photometry of M-dwarfs in SDSS is corrupted by a well-documented red leak problem, so we omit these filters.

Table 3. Summary of the Photometric and Astrometric Properties of the LSPM J1112+7626 System

Parameter	Value
α2000a,b	$11^{{\rm h}}12^{{\rm m}}42\mbox{$.\!\!^{\mathrm s}$}35$
δ2000a,b	+76°26'563
μαcos δb	0131 yr−1
μδb	−0108 yr−1
gc	15.618 ± 0.010
rc	14.195 ± 0.012
zc	12.027 ± 0.012
ICd	12.14 ± 0.05
J2MASSe	10.591 ± 0.023
H2MASSe	9.989 ± 0.021
K2MASSe	9.728 ± 0.017

Notes. ^aEquinox J2000.0, epoch 2000.0. ^bFrom Lépine & Shara (2005). ^cFrom SDSS Data Release 7. Note that these uncertainties do not account for variability and are therefore underestimates. ^dFrom Section 2.4. ^eWe quote the combined uncertainties from the 2MASS catalog, noting that the intrinsic variability of our target means in practice that these are underestimates.

Download table as:  ASCIITypeset image

The moderately high proper motion of our target allows the contribution of any background objects to the system light to be constrained using previous epochs of imaging. In Figure 2, we show imaging data with the position and size of the MEarth photometric aperture (which is representative of the aperture sizes used for all the photometry) overlaid. The DSS1 image indicates the aperture is unlikely to contain any significant flux from background objects, and the SDSS image (which has the highest angular resolution) does not appear elongated.

Zoom In Zoom Out Reset image size

A single epoch intermediate-resolution (R ≃ 1200) spectrum of our target was obtained using the FLWO 1.5 m Tillinghast reflector and FAST spectrograph (Fabricant et al. 1998) on UT 2011 March 16. The default setting of a 300 g mm⁻¹ grating and 3'' slit was used to obtain two exposures of 300 s. The spectra were reduced using standard procedures for long slit spectra in IRAF⁹ (Tody 1993), using an HeNeAr arc exposure taken in between the two target exposures to ensure accurate wavelength calibration. The two-dimensional spectra were divided by a flat field before extraction to remove the majority of the fringing seen at the reddest wavelengths. Relative flux calibration was performed using the spectrophotometric standard BD+26 2606. We show the resulting spectra in Figure 3.

Zoom In Zoom Out Reset image size

Using the TiO5 band index defined by Reid et al. (1995), we find a composite spectral type of M3.7 (note that the uncertainty in this value is at least 0.5 subclass), and using The Hammer (Covey et al. 2007), we find M3-M4 (visual inspection indicates the observed spectrum is closer to M4). The spectral type and the measured colors are in good agreement using the empirical colors of Galactic disk M-dwarfs from Leggett (1992) for the Johnson/Cousins/CIT system and West et al. (2005) for the SDSS/Two Micron All Sky Survey (2MASS).

It is also clear from Figure 3 that no Hα line is seen in the spectrum. Such behavior is typical of inactive field dwarfs at these spectral types when observed at low resolution, where the equivalent width of the Hα absorption line becomes close to zero at M4-M5 (e.g., Gizis et al. 2002; see their Figure 5, noting that the TiO5 index was 0.42 for LSPM J1112+7626). The presence of Hα emission in M-dwarf spectra is usually indicative of high activity levels.

The relative strengths of the TiO and CaH bands in optical spectra have been used to identify and classify M-subdwarfs (Gizis 1997), and more generally proposed as a method for estimating metallicities of M-dwarfs (e.g., Woolf & Wallerstein 2006). While detailed calibrations (e.g., using band indices) do not yet appear to be available, we performed a visual comparison between Figure 3 and the spectra shown in Gizis (1997). The strong TiO bands seen in our target indicate that it is not extremely metal poor. We therefore estimate it has metallicity [Fe/H] > −1 and is probably closer to solar (or greater).

High-resolution spectroscopic observations were obtained using the TRES fiber-fed échelle spectrograph on the FLWO 1.5 m Tillinghast reflector. We used the medium fiber (23 projected diameter), yielding a resolving power of R ≃ 44,000.

Observations commenced on UT 2010 May 1, and continued until UT 2011 April 13, with 43 epochs acquired. Total exposure times ranged from 45 to 80 minutes, broken into sequences of three or four individual exposures. The spectra were extracted using a custom built pipeline designed to provide precise radial velocities for échelle spectrographs. The procedures are described in more detail in Buchhave et al. (2010) and will be summarized briefly below.

In order to effectively remove cosmic rays, the sets of raw images were median combined. We used a flat field to trace the échelle orders and to correct the pixel-to-pixel variations in CCD response and fringing in the red-most orders, then extracted one-dimensional spectra using the "optimal extraction" algorithm of Hewett et al. (1985) (see also Horne 1986). The scattered light in the two-dimensional raw image was determined and removed by masking out the signal in the échelle orders and fitting the inter-order light with a two-dimensional polynomial.

Thorium–Argon (ThAr) calibration images were obtained through the science fiber before and after each stellar observation. The two calibration images were combined to form the basis for the fiducial wavelength calibration. TRES is not a vacuum spectrograph and is only temperature controlled to 0.1 °C. Consequently, the radial velocity errors are dominated by shifts due to pressure, humidity, and temperature variations. In order to successfully remove these large variations (>1.5 km s⁻¹), it is critical that the ThAr light travels along the same optical path as the stellar light and thus acts as an effective proxy to remove these variations. We have therefore chosen to sandwich the stellar exposure with two ThAr frames instead of using the simultaneous ThAr fiber, which may not exactly describe the induced shifts in the science fiber and can also lead to bleeding of ThAr light into the science spectrum from the strong argon lines, especially at redder wavelengths. The pairs of ThAr exposures were co-added to improve the signal-to-noise ratio, and centers of the ThAr lines found by fitting a Gaussian function to the line profiles. A two-dimensional fifth-order Legendre polynomial was used to describe the wavelength solution.

To model the system, we used the light curve generator from the popular jktebop code (Southworth et al. 2004a, 2004b), which is in turn based on the EB orbit program (ebop; Popper & Etzel 1981). This model approximates the stars using biaxial ellipsoids (Nelson & Davis 1972), which is appropriate for well-detached systems such as the present one. We modified the model to perform the computations in double precision, to account for light travel time, which is necessary due to the large semimajor axis, and to incorporate a simple prescription for the effect of photospheric spots (discussed in Section 4.3). The light curve model was coupled with a standard Keplerian orbital model for the radial velocities (also including the light travel time). We correct the orbital period from the heliocentric frame to the system barycenter (Griffin et al. 1976) when computing the physical parameters, following the usual convention of reporting the orbital period P as measured in the heliocentric frame (i.e., without the correction) in the tables.

We fit all of the Clay, Hankasalmi, and 2010–2011 season MEarth data, the radial velocities, and the spectroscopic light ratio simultaneously. We also include the single primary eclipse of UT 2010 April 3 to improve the constraint on the ephemeris, assuming the bandpass was approximately the same as the I_Bessell filter (this makes little difference in practice due to the sparse coverage). The remainder of the MEarth 2009–2010 season data were omitted as these are of quite poor quality and were taken in a different bandpass from the remainder, so do not usefully constrain the model after accounting for the extra free parameters which must be added in order to fit them.

Table 5 summarizes the assumptions made in our light curve model. The light curve parameter set generally follows jktebop and has been chosen to minimize correlations between parameters where possible. We add a term k_l(X − 1) to the model magnitudes, where X is airmass, and l is an index for the different light curves, to account for differential atmospheric extinction between the target and comparison stars, presumably due to color dependence. A different k_l value was allowed for each light curve l. This has only been used for fitting the Clay secondary eclipse curves and 2010–2011 MEarth data (see Section 2.2), where it seems to be necessary to include this term to obtain an acceptable fit to the out-of-eclipse portions. We are unsure as to the exact origin of the baseline curvature in the Clay data, and note that the choice of a linear function of airmass is arbitrary and may not be correct. This is further evidenced by the inconsistent sign and magnitude of this term in Table 6. We suggest the best approach to rectifying this problem would be to obtain higher quality secondary eclipse observations.

Gravity darkening coefficients were fixed at values appropriate for convective atmospheres (Lucy 1967), although gravity darkening is unimportant with the precision of the present light curves at such low rotation rates. The reflection effect was computed rather than fitting for it, although again the large semimajor axis means it is unimportant.

We used a square root limb darkening law (Diaz-Cordoves & Gimenez 1992) of the form

$$\begin{equation} {I_\lambda (\mu)\over {I_\lambda (1)}} = 1 - u (1 - \mu) - u^{\prime } (1 - \sqrt{\mu }), \end{equation} \tag{ 1 }$$

where I_λ(μ) is the specific intensity, and μ = cos θ, where θ is the angle between the surface normal and the line of sight. As discussed by van Hamme (1993), the square root law is a better approximation to the specific intensity distribution given by model atmospheres of late-type stars in the NIR than the other common two-parameter laws (quadratic and logarithmic) implemented in jktebop. We also verified this by comparing the Claret (2000) four-parameter law (which we assumed to be the best representation of the original phoenix model) with the two-parameter laws for typical M-dwarf stellar parameters.

For our baseline model, we fixed the limb darkening coefficients to values for the I_C filter derived from phoenix model atmospheres by Claret (2000), using T_eff = 3130 K for the primary, and T_eff = 3010 K for the secondary, log g = 5.0 and solar metallicity. We assume the same limb darkening law applies to all the I filters (the error introduced by this assumption is examined in Section 4.5).

No evidence was found for any background sources in the photometric aperture (see Section 2.5), and the radial velocity observations do not show additional lines or velocity drifts that might indicate a comoving, physically associated companion with the EB system. We therefore assume no third light (L₃ = 0), but see also Section 4.6.

To derive parameters and error estimates, we adopt a variant of the popular Markov chain Monte Carlo analysis frequently applied in cosmology and for analysis of exoplanet radial velocity and light curve data (e.g., Tegmark et al. 2004; Ford 2005). We use adaptive Metropolis (Haario et al. 2001) rather than the standard Metropolis–Hastings method with the Gibbs sampler. The adaptive method allows adjustment of the proposal distributions while the chain runs and is thus simpler and faster to operate, obviating the need for manual tuning of the proposal distributions typical of the more conventional methods. This method uses a multivariate Gaussian proposal distribution to perturb all parameters simultaneously.

Compared to the simpler Monte Carlo and bootstrapping methods (e.g., as implemented in jktebop), these methods have an important advantage for the present case where a heterogeneous set of light curves and radial velocities must be analyzed: it is possible to estimate the appropriate inflation of the observational uncertainties self-consistently and simultaneously with the other model parameters (e.g., Gregory 2005). This allows the correct relative weighting of the various observational data sets to be decided essentially by their residuals from the best-fitting model, and the uncertainty in this weighting to be propagated to the final parameters. Of course, this assumes the correct model has been chosen for the data, so these methods must be used carefully.

In the Monte Carlo simulations, we used the standard Metropolis–Hastings acceptance criterion, accepting the new point with probability:

$$$ P_{\rm acc} = \left\lbrace \begin{array}{@{}l} \dsty 1, P(M_{i+1}|D) \ge P(M_i|D) \\ \ms \dsty {P(M_{i+1}|D)\over {P(M_i|D)}}, P(M_{i+1}|D) < P(M_i|D), \end{array} \right. $$$

where P(M_i|D) and P(M_{i + 1}|D) are the posterior probabilities of the previous and current points, respectively, M_i represents the ith model, and D the data. The previous point in the chain was repeated if the new point was not accepted (this is required to correctly implement the method). The appropriate ratio of posterior probabilities ("Bayes factor") was

$$\begin{equation} {P(M_{i+1}|D)\over {P(M_i|D)}} = {P(M_{i+1})\over {P(M_i)}}\ \exp (-\Delta \chi ^2/2), \end{equation} \tag{ 2 }$$

where the factors in the right-hand side are the ratio of prior probabilities and the usual (least squares) likelihood ratio, respectively.

Priors assumed for each parameter are detailed in Table 5 and are chosen to be uninformative. In most cases, this choice is a uniform improper prior (labeled "uniform" in the table). For the orbital inclination, we adopt a uniform prior on cos i, which produces an isotropic distribution of orbit normals, and for the eccentricity, we use a uniform prior on e and ω, rather than on the jump parameters ecos ω and esin ω (the latter would produce a prior in e, ω proportional to e; e.g., Ford 2006).

As discussed by Gregory (2005), the appropriate choice of uninformative prior for "scale parameters" such as radial velocity amplitude or orbital period is a Jeffreys prior (or a modified Jeffreys prior for parameters which can be zero). This enforces scale invariance (equal probability per decade). In the present case, some of these parameters (particularly the period) are so well-determined the choice of prior is unimportant, so for simplicity we have only used non-uniform priors on the less well-determined parameters where the prior is important.

As discussed earlier, and in Gregory (2005), it is possible to fit for the appropriate inflation of the observational uncertainties simultaneously with the model parameters. This has been done via the s_k parameters shown in Table 5. It is conventional to apply this inflation by multiplication for light curves, and for radial velocity (RV) by adding an extra "jitter" contribution in quadrature. We follow these conventions here. Note that additional factors appear in the likelihood ratio of Equation (2) when doing this, which have been subsumed into the Δχ² in the equation.

For the RV, it is difficult to estimate reliable uncertainties from cross-correlation functions, so instead we derive them during the simulation. We weight each point by C²_peak (the square of the peak normalized cross-correlation value; see Table 4), where we find C_peak is approximately proportional to the signal-to-noise ratio in typical observations at similar signal to noise and peak correlation values to those seen in this work. This procedure is essentially equivalent to photon weighting, except the cross-correlation accounts for all sources of uncertainty rather than only photon noise. The uncertainty corresponding to a peak correlation of unity is derived during the simulation, and we allow separate values of this quantity for each star, labeled s₁ for the primary, and s₂ for the secondary, in the tables. The resulting per-data-point uncertainties used in the fitting are given by σ_ij = s_i/C_{peak, j} for the ith star and jth radial velocity point.

All Monte Carlo simulation results reproduced in this paper are derived from chains of 2 × 10⁶ points in length, where the first 10% of the points were used to initialize the parameter covariance matrix for the adaptive Metropolis method (starting from the initial parameters and covariance matrix derived using a Levenberg–Marquardt minimization¹⁰), and the next 40% were used to "burn in" the chain. This was found to be sufficient to ensure the chains were very well converged, and all of these first 50% of the points were discarded, leaving the remainder (10⁶ points) for parameter estimation. Correlation lengths in all parameters were <200 points.

We report the median and 68.3% central confidence intervals as the central value and uncertainty for all parameters. Reduced χ² values for all the model fits were unity by construction.

In LSPM J1112+7626, the presence of spots is clearly indicated by the observed out-of-eclipse modulation described in Section 2.4.

Spots present severe difficulty for the analysis of EB light curves because they cause systematic errors in the measurements derived from the eclipses, principally their depths, which determine the geometric parameters J, (R₁ + R₂)/a, and i (and by extension, R₂/R₁). Spots occulted during eclipse temporarily reduce the depth as they are crossed, and spots not occulted during eclipse increase it, because the surface brightness of the photosphere under the eclipse chord is greater than would be inferred from the out-of-eclipse baseline level.

The difficulty in solving for physical system parameters arises because the true spot distribution is not usually known. Out-of-eclipse modulations, and in non-synchronized systems, modulation of the eclipses themselves, are sensitive only to the longitudinal inhomogeneity in the spot distribution, except in special cases where spots are crossed during eclipse and cause detectable deviations from the usual light curve morphology (this can be difficult to detect in systems with radius ratios close to unity, because the deviations then last a large fraction of the eclipse duration). Any longitudinally homogeneous component, such as a polar spot viewed equator-on, or a homogeneous surface coverage of small spots, cannot usually be detected from light curves.

This can cause undetected systematic errors in the radii and effective temperatures derived from the light curves. Morales et al. (2010) discuss this issue in some detail for the best-measured literature EBs, finding this effect could produce up to 6% systematic errors in the radii for stars with spots of 30% filling factor, when the spots are concentrated at the pole.

The influence of the spectroscopic light ratio is not commonly discussed with regard to spots, but this is important because it provides complementary information to the photometry, on the difference in spot coverage between the two stars. For example, consider the case of a large coverage of longitudinally homogeneous, polar spots that are not crossed during eclipse. If these spots are distributed on both stars, the light ratio is unaffected but both eclipses will appear deeper. However, if the spots are on only one star, the effect on the photometry is the same, but the light ratio is altered because the spotted star appears darker. This will produce a discrepancy between the light ratio derived only from the light curve parameters, and the spectroscopic value. More generally, it is possible to place constraints on the difference in the overall spot coverage between the binary components, using the spectroscopic information.

One unusual feature of the present system among low-mass EBs showing signs of spots is the spin and orbit appear to not be synchronized. Because of this, it is possible to measure eclipses at different rotational phases, and thus different spottedness of the visible stellar hemispheres. While this does not necessarily resolve the difficulty of determining the longitudinally homogeneous component of the spot distribution, it does provide additional information on the longitudinally inhomogeneous component, specifically which star the spots are located on. This information is not usually available if only out-of-eclipse modulations are seen. It is common to assume the spots are on both binary components, which may be reasonable for near-equal mass systems at short periods where tidal synchronization is expected to have occurred, but this is not a reasonable assumption for LSPM J1112+7626, where only one out-of-eclipse modulation is seen. The observed rotational evolution of M-dwarfs (e.g., Irwin et al. 2011) indicates that it is extremely unlikely the two binary components could have the same rotation period (and phase) by chance in the absence of tidal effects, which are not expected.

Unfortunately, at the time of writing, only a single primary eclipse is available, so we are unable to take full advantage of the non-synchronized spin and orbit at present. Also, while multiple secondary eclipses were observed, and residuals from our best-fitting model (assuming no changes in eclipse shape or depth) are seen, it is not clear if many of these are simply the result of systematic errors in the photometry, given that equally large deviations are seen out of eclipse. Constraining the spot properties by this method will be an important area for future work, and observing in multiple bandpasses would be advantageous as it may provide information on the spot temperatures.

4.3.1. Spot Model

Conventionally, a Roche model such as the one implemented in the popular Wilson–Devinney program (Wilson & Devinney 1971) would be used to model a system with spots as these usually include a circular spot model and perform the necessary surface integrals over the stars; however, the treatment of proximity effects is completely unnecessary in the present system, and these models are extremely computationally intensive, especially for systems with eccentric orbits, which would make a detailed Monte Carlo based error analysis prohibitive.

It is also not clear that the circular spot model with a small number of spots is realistic for late-type dwarfs. Observed light curves very rarely show the characteristic "eclipse-like" features with flat baseline at maximum light, as would be predicted for a small number of near-equatorial spots. This indicates either that these objects have very large, polar spots (e.g., Rodono et al. 1986), such that some of the spot is always in view to the observer to produce the continuous photometric modulations seen in light curves, or that the surfaces of these objects have many small spots (e.g., Barnes & Collier Cameron 2001; Barnes et al. 2004) with the photometric modulations arising from a longitudinal inhomogeneity of the spot distribution. Modeling an ensemble of spots, as in the latter case, in detail would be prohibitive, as even single spot models are usually degenerate given limited light curve data. In our case, the degeneracies in fitting spot models would be further exacerbated by the availability of only single-band light curve information, meaning spot temperature and size would be essentially degenerate.

Given these difficulties, we take an alternative approach for the solution presented in this paper. Rather than trying to model the unknown spot distribution in detail, we simply attempt to incorporate the effect of spots into the uncertainties on the final parameters. We consider two models: (1) spots on the non-eclipsed part of the photosphere on both stars and (2), a case intended to approximate the effect of spots covering the whole photosphere, again for both stars. In both cases, the models are required to reproduce the observed out-of-eclipse modulation. It is necessary to run two models in each spotted case because it is not known which star the out-of-eclipse modulation originates from.

We assume a simple sinusoidal form for the modulations, which appears to be an adequate description of the available light curve data. The functional form adopted was

$$\begin{equation} {\Delta L_i\over {L_i}} = a_i \sin \left({2 \pi F_i t\over {P}}\right) + b_i \cos \left({2 \pi F_i t\over {P}}\right) - \sqrt{a_i^2+b_i^2}, \end{equation} \tag{ 3 }$$

where t is time, L_i is the light from star i, a_i and b_i are constants expressing the amplitude and phase of the modulation, and F_i is the "rotation parameter," the ratio of the rotation frequency to the orbital frequency. This expression is constructed to yield ΔL_i = 0 at maximum light, corresponding to the conventional assumption that no spots are visible when the star is brightest, which gives the minimum spot coverage necessary to reproduce the observed modulations.

jktebop computes the final system light (and thus the change in magnitude) by summing the out-of-eclipse light and then subtracting the eclipsed light. By modulating only the out-of-eclipse light in this summation using Equation (3), hypothesis (1) can be implemented and hypothesis (2) can be implemented by also modulating the eclipsed light (changing the surface brightness under the eclipse chord) when the eclipsed star is star i.

4.3.2. Results

The effect of applying our spot model on the eclipses is shown in Figures 4 and 5, which were calculated using the parameters of LSPM J1112+7626, comparing a model with spots to a model with identical physical parameters, but without spots. This figure demonstrates that observing multiple eclipses at high precision would allow the star hosting the spots to be identified.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We fit all four spotted models to the full data set for LSPM J1112+7626 using the procedures already described. The results are given in Table 6, and Figures 6–8 show the data with a representative model, corresponding to hypothesis (1), the "non-eclipsed spots" model, on the primary, overplotted.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The 2010 November 19 partial secondary eclipse from MEarth (numbered 5 in the figure) provides the most discriminating power between the various models, particularly taken in conjunction with the 2011 May 2 event (numbered 9), which dominates the fit. As shown in Figure 7, these were taken at opposite extremes of the out-of-eclipse variation, with the 2010 November 19 event at maximum light and the 2011 May 2 event close to minimum light. Examining the sizes of the s parameter in Table 6, the "eclipsed spots" model on the secondary, which is the only model presenting a significant detectable deviation in the secondary eclipses (as shown in Figure 5) is disfavored by the data (see Figure 9), having an s_{MEarth-20101119} value larger by almost 3σ than the "non-eclipsed spots" model on the secondary. Indeed, the observations indicate that this eclipse needs to be made as shallow as possible in the model in order to produce a good fit, possibly even slightly shallower than the "non-eclipsed spots" model is able to produce (the discrepancy may arise from the imperfect modeling of the out-of-eclipse variation by a sinusoid). The depth of this eclipse therefore argues that if the spots causing the out-of-eclipse modulation are on the secondary, they should be at latitudes not crossed by the eclipse chord.

Zoom In Zoom Out Reset image size

It is interesting to note that, because the radius ratio is quite close to unity, the range of allowed combinations of spot latitude and inclination angle for the secondary that would explain the observed out-of-eclipse modulations is therefore quite narrow, compared to the much wider range of possibilities on the primary that are still allowed. This mildly argues in favor of spots on the primary as being the more likely explanation, purely from geometry.

Nevertheless, with the presently available data, the other models do not show any clear disagreement with the observations, so we consider (from a purely observational point of view) that all three are equally likely at present.

4.3.3. Effect of Additional Spots

A concern, when combining light curves taken with different instruments, is the effect of any difference in the photometric bandpasses. In the present case, where the primary and secondary eclipses were (by necessity) measured using different instruments, this predominantly affects the ratio of the eclipse depths, and thus the parameter J, the ratio of central surface brightnesses.

In order to determine the approximate size of the bandpass mismatch, we measured the difference in magnitude between LSPM J1112+7626 and a nearby, much bluer comparison star, at $11^{{\rm h}}12^{{\rm m}}12\mbox{$.\!\!^{\mathrm s}$}15$+76°27'338 (position from 2MASS; this is the bright star seen near the upper left corner of Figure 2). The latter star has (J − K)_2MASS = 0.27 and is thus probably in the F spectral class. It is non-variable at the precision of the MEarth data (rms scatter 0.003 mag).

Our measured magnitude differences were ΔI(Hankasalmi) = 0.36, ΔI(Clay) = 0.35, and ΔI(MEarth) = 0.55. The excellent agreement between the two I_Bessell filters (especially in light of the observed out-of-eclipse variations of our target) indicates they are most likely an extremely good match.

Assuming a simple linear scaling, we estimate the change in J corresponding to the mismatch between the two I_Bessell filters is approximately 1/10 of that between the I_Bessell and I_MEarth filters, or δJ ≈ 0.003. This is unimportant compared to the other sources of uncertainty (see Table 6).

4.5.1. Differences between I_Bessell and I_C

4.5.2. Errors in the Atmosphere Models

We now examine our assumption of no third light. While the proper motion evidence (see Figure 2) argues it is unlikely there are any background stars contributing at a significant level in the photometric aperture, the relatively low proper motion of our target does not allow this possibility to be completely eliminated given the limited angular resolution of the first epoch imaging data. Common proper motion companions are also permitted at intermediate separations below the angular resolution of the SDSS data, but still in wide enough orbits or with mass ratios much less than unity, where they would not be detected by the presence of additional lines in the spectra or radial velocity drift over the approximately one year baseline available.

Therefore, to examine any constraints on third light which can be placed directly from the light curves, and the effect on the derived parameters, we ran an additional set of Monte Carlo simulations using the basic "non-eclipsed spots" model on the primary discussed in Section 4.3.2 and Table 6, allowing L₃ to vary. For simplicity, a uniform prior was assumed. We find L₃ < 0.029 at 95% confidence. As discussed by Nelson & Davis (1972), the dominant parameter mimic is between third light and inclination, and our results confirm this, finding the uncertainties in cos i and the central surface brightness ratio parameters were slightly inflated. The other light curve parameters and radii derived from these chains are indistinguishable within the uncertainties from the results where L₃ was not varied.

Of the sources of uncertainty we have considered, spots dominate. Since the stars hosting the spots and the spot configuration are mostly unknown, we adopt the union of the six spotted solutions discussed in Section 4.3.3, giving them equal weight. The final parameter estimates and their uncertainties are reported in Table 10, derived from the posterior samples produced by the Monte Carlo simulations, adopting the median and 68.3% central confidence intervals as the central value and uncertainty (see Section 4.2).

We compute UVW space motions for LSPM J1112+7626 using the method of Johnson & Soderblom (1987), with the distance from the EB solution, γ velocity from Table 6, and position and proper motions from Table 3. We find this object is in the old Galactic disk population, following the method and definitions in Leggett (1992).

For very precise work, even matters as seemingly trivial as physical constants can be important, so we briefly discuss our assumptions in this regard. We adopt the 2009 IAU values of G M_☉ and the astronomical unit (see the Astronomical Almanac 2011), and follow Cox (2000) in adopting the value of the solar photospheric radius from Brown & Christensen-Dalsgaard (1998).¹¹ For the solar effective temperature, we use the value of 5781 K from Bessell et al. (1998), which is based on measurements of the solar constant by Duncan et al. (1982).

The final constant, the solar bolometric absolute magnitude, is a more thorny issue, as discussed by Bessell et al. (1998), and extensively by Torres (2010). When using tables of bolometric corrections, it is vital to adopt the consistent value of this quantity in accordance with the table. We therefore use the appropriate values in each case in the following section.

In order to estimate effective temperatures (and thus, bolometric luminosities, and the distance), an external estimate of the effective temperature of one of the binary components is required, in addition to bolometric corrections. For M-dwarfs, there are large systematic uncertainties in the effective temperature scale (e.g., Allard et al. 1997; Luhman & Rieke 1998), with disagreement at the few hundred degrees Kelvin level among different authors. Many of the early difficulties were due to the lack of model atmospheres (which are usually needed to integrate the full spectral energy distribution from the available measurements to obtain L_bol), and significant improvement on this front was made in the 1990s. The availability of a much larger sample of angular diameter measurements from interferometry should further improve the situation as this provides a much more direct method to estimate T_eff for single stars, but there are relatively few temperature scales available at the present time using this information.

We show results from two different inversions to illustrate the typical range of parameters. Both scales we use were derived with model atmospheres rather than blackbodies, and in both cases we use the measured I_C − J_2MASS, I_C − H_2MASS, and I_C − K_2MASS colors from Table 3 in conjunction with the component radii and I_C-band light ratio from the adopted EB solution. The results from these different pairings of bandpasses were found to be consistent within the uncertainties, so we report the union of the results.

The two sets of effective temperature and bolometric corrections adopted were from Casagrande et al. (2008), which is a recent determination using interferometric angular diameter measurements to derive T_eff, and from Leggett et al. (1996), which is the basis for several recent works and tabulations of the effective temperature scale, including Luhman & Rieke (1998), Luhman (1999), and Kraus & Hillenbrand (2007).

For the Leggett et al. (1996) scale, we fit polynomials to their tabulated data, omitting three objects these authors found to be metal poor from our fits: Gl 129, LHS 343, and LHS 377. We follow Leggett et al. (1996) in adopting a systematic uncertainty of 150 K, which appears to be consistent with the differences we find between the two determinations.

Our results for LSPM J1112+7626 are reported in Table 10.

We first limit our comparisons to the mass–radius plane. If the metallicity is solar, the sum of radii is inflated at the 7σ level compared to the models, even for 10 Gyr age, which produces the largest model radii. The discrepancies are not as significant in the individual radii, with the primary 2σ and the secondary 3σ larger than the model predictions. We also note that the "eclipsed spots" model on the secondary, which we discarded in Section 4.3, does not solve the radius discrepancy, having a secondary radius (and sum of radii) even further from the model predictions.

Youth is unlikely as an explanation given the old disk kinematics of our target. We estimate an age of approximately 120 Myr at solar metallicity would be needed to bring the radius sum into agreement with the models, yet such objects should also rotate extremely rapidly and exhibit strong Hα emission (e.g., Terndrup et al. 2000; Hartman et al. 2010, and references therein), neither of which are seen.

Metallicity may provide a viable explanation for the observed radii given the present uncertainties (both in the observations and the models, where for the latter it is possible the effect of metallicity on radius is underpredicted; see Section 1 and López-Morales 2007), although the object would most likely need to be extremely metal rich. This would place it in the tail of the metallicity distribution for the old Galactic disk population, which is likely to peak around [M/H] ≈ −0.5 (e.g., Leggett 1992, and references therein).

Given that metallicity has difficulty reproducing radii in other well-observed systems (see Section 1) and the kinematic evidence, we proceed to examine other possible causes of the inflated radii.

It is tempting to suggest that the inflation in LSPM J1112+7626 results from elevated activity levels in one or both components. While the lack of Hα emission, unless this is rendered undetectable by the limited signal-to-noise ratio of our spectra, may be a difficulty for this scenario to explain, it should be noted that more modest levels of activity can result in weakened Hα emission, or even absorption, but may still produce sufficient spot coverage to explain the inflation.

This system straddles (in mass) the full convection limit, where a substantial increase in the observed activity lifetimes for M-dwarfs occurs (e.g., West et al. 2008). This would argue that the secondary is more likely to be active and thus the source of the inflation. The secondary lines in the spectrum are also more difficult to observe due to its lower luminosity, which would make an Hα line easier to hide in the noise. Identification of the star hosting the out-of-eclipse modulations may shed additional light on which of the components might be responsible for the inflation.

As discussed by Morales et al. (2010), the radius inflation inferred from EB analyses often implies extremely high spot filling factors. We estimate that ρ = 1.05 corresponds to β = 0.1. Assuming a spot temperature contrast of T_{eff, spot}/T_{eff, phot} = 0.9, this would require a filling factor of approximately 30%. This is much larger than needed to produce the observed out-of-eclipse modulation and would mean there is a substantial longitudinally symmetric spot component, if the inflation is indeed caused by spots. We showed in Section 4.3.3 that such a level of spots does slightly modify the component radii if they are not eclipsed, reducing R₂ and increasing R₁. In practice, this may slightly reduce the required filling factor depending on where the spots are located, but unless the spot coverage is much larger than we assumed, the systematic errors in the radii resulting from such spots are still insufficient to explain the observations without invoking inflation.

Finally, on suggestion of the referee, we note that this object may indicate the need to revisit the equation of state for very low mass stars, which was discussed as a possibility to explain the radius inflation before the activity hypothesis gained prevalence (e.g., Lopez-Morales 2004; Torres & Ribas 2002). The equation of state was a significant source of difficulty in early work attempting to model very low mass stars (Chabrier & Baraffe 1995), and although substantial progress on this front has been made (see the reviews by Chabrier & Baraffe 2000; Chabrier et al. 2005), it is still an open question.

Comparisons in only the mass–radius plane ignore important information contained in the effective temperatures. The correct physical model must explain all of the observations simultaneously, so we now proceed to examine the temperatures. This comparison is more problematic than in mass–radius for a number of reasons. The main observational issues are the difficulty of constraining effective temperatures (see Section 4.8) and metallicities for the EBs.

As noted in Section 1, unlike the radius, the effective temperature predicted from models does depend quite strongly on metallicity, because the bolometric luminosity is a function of metallicity. As the metallicity is decreased, the model bolometric luminosity and effective temperature increase, and the radius decreases by a small amount. While metallicity does complicate the interpretation of the mass–T_eff diagram, this is a key argument against metallicity as the explanation for all of the inflated radii in the EB sample, as discussed already in Section 1. In this regard, CM Dra is puzzling as several authors have claimed this object is metal poor, which would further exacerbate the effective temperature discrepancy with the models.

Given the present uncertainties, it is not clear if the effective temperature difference between LSPM J1112+7626 and the models is significant. However, it does lie at the low-temperature end of the EB results, and this lends support to the hypothesis that this object is more metal rich than the other EBs.

It is clear that LSPM J1112+7626 offers one of the best prospects to observationally test the causes of inflated radii in EBs, but further observations are needed. The most fruitful avenues would be to pursue a determination of the metallicity, improvement of the system parameters with a particular focus to investigating potential systematics in the radii (e.g., due to spots), and constraining the activity levels in the components by independent means (for example, X-ray emission, or the Caii H and K lines).

While we have mentioned the non-synchronized spin and orbit, we were unable to take full advantage of this property with the available observational data. As shown in Section 4.3, the effect of spots should be larger on the primary eclipse, so we advocate an intensive campaign on these events, preferably to obtain multiple, complete eclipses with the same telescope and detector system at a range of rotational phases and in multiple wavelengths. The primary eclipses are observable from Europe and Scandinavia during the winter 2011–2012 season, where excellent observational facilities are available. We remind observers that good-quality out-of-eclipse monitoring will be important for the interpretation of the eclipses themselves, and the spot distributions are likely to evolve (indeed, Figure 7 indicates that we may have already seen this in the existing data), so the monitoring must be contemporaneous, and preferably also taken at multiple wavelengths to constrain the spot temperatures.