GPI Spectroscopy of the Mass, Age, and Metallicity Benchmark Brown Dwarf HD 4747 B

Late-type dwarfs are assigned a spectral type by comparison of near-infrared spectra to that of spectral standards (Burgasser et al. 2006a). We employed this method in Crepp et al. (2015), and use it again here for the spectrum of HD 4747 B. We compare the spectrum of HD 4747 B to 308 L and T dwarf spectra from the SpeX prism library, supplemented by spectra from Filippazzo et al. (2015). Most spectra are from Burgasser et al. (2004, 2006a, 2008, 2010), Burgasser (2007), Chiu et al. (2006), Cruz et al. (2004), Kirkpatrick et al. (2011), Liebert & Burgasser (2007), Looper et al. (2007), Mace et al. (2013), Mainzer et al. (2011), and Sheppard & Cushing (2009).

We begin by binning and trimming comparison spectra to match the resolution and wavelength coverage of the HD 4747 B GPI spectrum. We compare the L and T dwarf spectra to that of HD 4747 B with a normalization constant and χ²-like goodness of fit statistic from Cushing et al. (2008). The χ² value is calculated for each binned L and T dwarf spectrum, incorporating errors from spectral templates and the spectrum of HD 4747 B. We perform this fit to the H- and K₁-bands separately and also to the "full" stitched HK₁ spectrum.

The top panels of Figures 2 and 3 display χ² values versus spectral type. The minimum χ² occurs in the early T dwarf range for the H-band-only spectrum. Meanwhile, the K₁-band and combined HK₁ comparisons both reach a minimum with an overall trend tending toward the L/T transition. Specific best matching stars are shown in Table 3. From this analysis, we qualitatively adopt a spectral type of T1 ± 2 noting the large scatter in the goodness-of-fit statistic between filter fits and fact that brown dwarfs with similar spectral features can have different physical properties anyhow.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Brown dwarf spectral types generally correlate with temperature. However, in order to infer physical properties, it is necessary to use methods complementary to the NIR spectral typing analysis performed above. We use the Saumon & Marley (2008) atmosphere model grid of synthetic spectra to fit HD 4747 B and infer its temperature. The models are designed to reproduce the L/T transition covering the temperature range needed for this application (Saumon & Marley 2008). They also incorporate an adjustable sedimentation factor (f_sed) and vertical mixing parameter (k_zz) to describe cloud clearing that nominally occurs during this transition. The models range from T = 900–1600 K in increments of ΔT = 100 K evaluated at surface gravities log(g) = 4.0, 4.5, 5.0, 5.5 (cgs units), sedimentation factor of f_sed = 0 (cloudless) or f_sed = 2.0 (cloudy), and eddy diffusion coefficient k_zz of either 0 (chemical equilibrium) or 4.0 (not in chemical equilibrium).

The fitting routine begins similarly to the spectral typing analysis above by binning and trimming the models to match the wavelength range and resolution of the spectrum of HD 4747 B. The model spectra are normalized to the companion spectrum via a normalization constant (Cushing et al. 2008). A χ² value is calculated for every binned model compared to the companion spectrum to determine a starting point on the model grid that results in a small burn-in phase for the next step in the Bayesian procedure. Probability distributions are generated using a Markov Chain Monte Carlo (MCMC) approach (Foreman-Mackey et al. 2013). The MCMC routine begins at the minimum χ² value parameter space point and linearly interpolates between model spectra flux.

We perform the fitting routine for the individual H- and K₁-bands separately as well as the combined HK₁ spectrum (Figure 3). The H- and K-band best-fit models were normalized to the individual band; the full spectrum was multiplied by the same normalization constant and overplotted on the companion spectrum for comparison. We also plot the combined H best-fit model at the same T_eff, log(g), and k_zz parameters with an f_sed of 0 (cloudless).

Zoom In Zoom Out Reset image size

We find that cloudless models are able to adequately fit only one filter data set at a time. Cloudy models do a much better job replicating the spectrum of HD 4747 B across the H- and K₁-bands, but discrepancies in overall infrared color remain (Figure 4). Similar results have been obtained for other brown dwarfs that show comparable spectra but distinct offsets in integrated flux (color; Leggett et al. 2002; Konopacky et al. 2016). This feature is commonly cited as a "secondary parameter" problem that involves missing physics in the models related to metallicity, gravity, cloud, and dust properties, unresolved binaries, the influence that electron conduction has on cooling, and other items that models may inherently assume, such as spherical geometry, hydrostatic equilibrium, nonmagnetic, nonrotating, etc. (Chabrier et al. 2000; Saumon & Marley 2008).

In any case, model fitting of our GPI spectrum suggests that HD 4747 B is a cloudy brown dwarf. Table 4 summarizes the results. Correlation matrices between parameters are shown in Figure 5. HD 4747 B shows an effective temperature that is consistent for each wavelength range to within ΔT ≈ 100 K. From this analysis, we adopt an effective temperature of T_eff = 1450 ± 50 K. This temperature is consistent with an L/T spectral type (Stephens et al. 2009; Filippazzo et al. 2015).

Zoom In Zoom Out Reset image size

Table 4.  HD 4747 B Model Fitting Results for H-band-only, K₁-band-only, and the Combined HK₁ Bands Using the Saumon & Marley (2008) Grids

Parameter	H	K1	HK1
Teff (K)	${1509}_{-141}^{+159}$	${1393}_{-150}^{+138}$	${1407}_{-140}^{+134}$
log(g) (cm s−2)	${5.2}_{-0.5}^{+0.5}$	${4.5}_{-0.5}^{+0.4}$	${5.2}_{-0.6}^{+0.5}$
fsed	0	0	${1.91}_{-0.18}^{+0.18}$
kzz	0	0	0

Download table as:  ASCIITypeset image

Model fits yield only relatively weak constraints on the companion surface gravity. The K₁-band fit suggests the lowest surface gravity ($\mathrm{log}g={4.5}_{-0.5}^{+0.4}\,\mathrm{dex}$), whereas the H-band fit ($\mathrm{log}g\,={5.2}_{-0.5}^{+0.5}\,\mathrm{dex}$) and full HK₁ spectrum ($\mathrm{log}g={5.2}_{-0.6}^{+0.5}\,\mathrm{dex}$) prefer a higher gravity (Table 4). In these models, k_zz = 0, while f_sed = 0 for the individual bands and f_sed = 1.9 for the joint HK₁ fit.

Since surface gravity is not always well constrained via model fitting of low-resolution spectra, we also infer $\mathrm{log}(g)$ from the companion dynamical mass (see Section 4.6) and radius estimated from evolutionary models. We use the independent gyrochronological age of the system, $\tau \,={3.3}_{-1.9}^{+2.3}\,\mathrm{Gyr}$, which was calculated in Crepp et al. (2016), based on the host star B − V value and estimated rotation period; this value includes scatter in the activity-rotation empirical relations from Mamajek & Hillenbrand (2008). We next use the Saumon & Marley (2008) models to determine an age-dependent radius, finding $R={0.089}_{-0.005}^{+0.006}\,{R}_{\odot }$. The radius in turn yields an estimate for the surface gravity, $\mathrm{log}(g)={5.33}_{-0.08}^{+0.03}$, which is based on the dynamical mass but is still model-dependent. Nevertheless, it offers a self-consistency check when comparing surface gravity with effective temperature, isochrones, and other parameters as is done in Section 4.7. For comparison, we find $R={0.087}_{-0.007}^{+0.005}\,{R}_{\odot }$ and $\mathrm{log}(g)={5.35}_{-0.08}^{+0.10}$ using the Baraffe et al. (2003) models. These surface gravity values are consistent with the spectral fitting results and also expected older age of the system, e.g., compared to young stellar clusters, moving groups, or associations.

Methane is a volatile compound whose presence as an absorption feature may be used to identify and classify brown dwarfs (Burgasser et al. 1999). Methane absorption in the H-band influences the near-infrared color of substellar objects and can help to distinguish between L dwarfs and T dwarfs (Geballe et al. 2002). To further assess our classification scheme, we analyze the H-band spectrum of HD 4747 B to search for methane by calculating a methane index, ${S}_{{\mathrm{CH}}_{4}}$,

$$\begin{eqnarray}&&{S}_{{\mathrm{CH}}_{4}}=\displaystyle \frac{{\int }_{1.560}^{1.600}{f}_{\lambda }d\lambda }{{\int }_{1.635}^{1.675}{f}_{\lambda }d\lambda },\end{eqnarray} \tag{ 1 }$$

where f_λ and dλ values are derived from the GPI spectrum. This technique was also used for an H-band methane analysis of the benchmark brown dwarf HD 19467 B with the Project 1640 IFS (Crepp et al. 2015). Brown dwarfs generally have methane indices ${S}_{{\mathrm{CH}}_{4}}\gt 1$. We find that HD 4747 B has ${S}_{{\mathrm{CH}}_{4}}=1.1\pm 0.1$. Using values for the spectral classification of brown dwarfs from Geballe et al. (2002), we infer a spectral type of L9.5-T2 (Figure 6). This range is consistent with our model fitting analysis and further supports the interpretation that HD 4747 B is an L/T transition object.

Zoom In Zoom Out Reset image size

At present, "atmospheric models remain largely untested at non-solar metallicities" (West et al. 2011). This is due in large part to a lack of available benchmark objects, HD 114762 B representing the latest available spectral type (d/sdM9) at [Fe/H] = −0.7 (Bowler et al. 2009). One of the goals of the TRENDS program is to detect benchmark companions with disparate metallicities. HD 4747 has a precisely measured metallicity of [Fe/H] = −0.22 ± 0.04 dex. However, for all intents and purposes, this value is very nearly solar, considering that metal content is expected to influence the emergent spectrum of brown dwarfs as a second-order effect. In particular, only cloudless models are publicly available that explore metallicity. Skemer et al. (2016) sample [M/H] at 0.0, 0.5, and 1.0 dex to study GJ 504 b, but custom models have not yet been developed for the L/T transition to our knowledge, so far only for transiting exoplanets (Morley et al. 2017).

As a provisional test of the influence of metallicity on HD 4747 B, we consider the recently developed Sonora model grids (Marley et al. 2017). Noting as a caveat that only cloudless models are currently available, we investigate qualitatively the spectral morphology and quantify the percentage difference between models sampled at solar metallicity and [Fe/H] = −0.5 dex. In the H-band, we find that the subsolar model offers a slightly better match to the GPI data though comparable to uncertainties. In the K₁-band, the shift in peak flux appears at shorter wavelengths, which arguably better reflects the spectral morphology of HD 4747 B, modulo the systematic flux offset noted previously. At most, the models differ by ≈15% considering that they straddle the expected metallicity of HD 4747 B, presuming it is comprised of the same material as its parent star. In comparison, the cloudy and cloudless models presented in Figure 4 differ by 41% on average across the H- and K-bands, indicating that the expected effect of subsolar metallicity on the spectrum of HD 4747 B is much less than that of its cloudy atmosphere. Theoretical work on the impact of metallicity on cloudy substellar objects is on-going and beyond the scope of this paper (M. Marley et al. 2018, in preparation).

We measure the position of HD 4747 B relative to its parent star to constrain the system orbit and dynamical mass. Observed in coronagraphic mode, the primary star remains hidden behind the occulting spot so its position must be inferred using the off-axis satellite spots described in Section 3. Following speckle suppression, we measure the centroid of the companion (and each injected fake companions) in each wavelength channel to determine the angular separation, position angle, and their uncertainties. We use the plate-scale of θ = 14.37 ± 0.34 mas pix⁻¹ from the GPI instrument calibration webpage based on data obtained by Konopacky et al. (2014) for Θ Ori using the orbital information from Magellan AO (Close et al. 2013).

Table 5 shows astrometric results for H and K₁ data respectively. Uncertainty in the plate scale is folded into the analysis by adding its contribution in quadrature with other uncertainties. We find that uncertainty in the position of HD 4747 B is dominated by speckles and residual atmospheric dispersion. Results from each filter are consistent with one another to within 1σ in angular separation and 2σ is position angle. The position of the companion with GPI are consistent with that of NIRC2 from the discovery article, which we also show for reference (Crepp et al. 2016). The companion appears to be moving in a counter-clockwise direction with an angular separation that is decreasing with time. Given that the GPI data was obtained only three months after the most recent NIRC2 data, we do not update the companion dynamical mass or orbit, which has a period of approximately 38 years. Performing a statistical (MCMC) analysis of the system astrometry, we find that posterior distributions yield comparable uncertainties as the original discovery article with the new GPI measurements.

Table 5.  Relative Astrometric Measurements of HD 4747 B Including Initial Discovery Data from NIRC2 (First Three Rows) and Two New GPI Observations (Bottom Two Rows)

JD-2,450,000	Instrument	Filter	ρ (mas)	Position Angle (°)	Proj. Sep. (au)
6942.8	NIRC2	Ks	606.5 ± 7.0	18004 ± 062	11.33 ± 0.17
7031.7	NIRC2	Ks	606.6 ± 6.4	18052 ± 058	11.34 ± 0.16
7289.9	NIRC2	L'	604 ± 7	1849 ± 09	11.3 ± 0.2
7380.6	GPI	K1	585 ± 14	1852 ± 03	10.93 ± 0.29
7381.5	GPI	H	583 ± 14	1844 ± 03	10.90 ± 0.29

Download table as:  ASCIITypeset image

We update the dynamical analysis presented in the HD 4747 B discovery paper (Crepp et al. 2016). We use the same methodology by jointly fitting the 18 years of stellar RV data with imaging astrometry from Table 5. GPI measurements extend the astrometric time baseline by ∼92 days. This is small compared to the expected P = 38 year orbit of the companion but, nevertheless, allows us to assess any systematic errors and perform consistency checks. Measurement uncertainties are comparable to the NIRC2 data set.

The python-based "MCMC Hammer" routine is used within a Bayesian framework to calculate posterior distributions for orbital parameters (K, P, e, ω, Ω, i, t_p), where K is the RV semi-amplitude, P represents the orbital period, e is eccentricity, ω is the longitude of periastron, Ω is the longitude of the ascending node, i is inclination, and t_p is time of periastron passage (Foreman-Mackey et al. 2013). Figure 7 shows results for the orbital parameter distributions, including nuissance parameters such as the global RV offset, differential RV offset (resulting from two distinct HIRES instrument settings), and RV jitter value. Parameter confidence intervals from the model are listed in Table 6.

Zoom In Zoom Out Reset image size

Table 6.  Orbital Parameters and Dynamical Mass of HD 4747 B

Orbital Fitting Results
Parameter	Value (68%)
K (m s−1)	${755.1}_{-11.6}^{+12.8}$
P (years)	${37.85}_{-0.78}^{+0.87}$
e	${0.740}_{-0.002}^{+0.002}$
ω (°)	${269.3}_{-0.6}^{+0.6}$
tp (year)	${1997.04}_{-0.02}^{+0.02}$
i (°)	${58.5}_{-5.2}^{+4.0}$
Ω (°)	${186.2}_{-1.7}^{+2.3}$
a [AU]	${12.3}_{-1.5}^{+1.5}$
Global offset (m s−1)	$-{214.4}_{-11.9}^{+10.9}$
Diff. offset (m s−1)	${13.2}_{-2.7}^{+2.8}$
Jitter (m s−1)	${4.0}_{-0.7}^{+0.9}$
mass [MJup]	${65.3}_{-3.3}^{+4.4}$

Download table as:  ASCIITypeset image

We find that most parameters remain very similar to the original analysis with the exception of the orbital inclination and semimajor axis, which are now both smaller and with reduced uncertainties. Although it is still early to confidently extract orbital information from astrometry data given the limited baseline, the result of a lower inclination and semimajor axis was actually predicted by Crepp et al. (2016) based on model-dependent estimates of the companion mass on the grounds of self-consistency. This result suggests that the GPI and NIRC2 measurements are compatible with one another at the ∼10 mas level. Figure 8 shows the most recent sky-projected orbit trajectories consistent with both imaging and RV data sets with predictions for future epochs.

Zoom In Zoom Out Reset image size

Finally, we use values directly from the MCMC chain parameters for K, i, e, and P to solve for the companion dynamical mass. To accommodate for uncertainty in the mass of the host star, we draw from a normal distribution centered on the spectroscopically determined mass from the discovery paper, ${M}_{* }=0.82\pm 0.04\,{M}_{\odot }$ (Crepp et al. 2016).¹² We numerically calculate the companion mass from each MCMC iteration to generate a posterior distribution. While the star mass may contain small, model-dependent systematic errors, this method ensures that uncertainty in the companion mass scales as $m\sim {M}_{* }^{2/3}$ making this contribution comparable to uncertainties that result from minimal astrometry orbital coverage.

Figure 9 shows results for the companion mass posterior before and after including GPI measurements. We find that the most likely (mode) companion mass has increased slightly since the discovery data, from to 60.2 ± 3.3 M_Jup to ${65.3}_{-3.3}^{+4.4}\,{M}_{\mathrm{Jup}}$, a consequence of the smaller inclination. The posterior distribution is still asymmetric and actually slightly broader than before, indicating that additional measurements are needed to further constrain the companion mass. Nevertheless, this information is sufficient to estimate the companion surface gravity, as was done in Section 4.2, and compare to evolutionary models.

Zoom In Zoom Out Reset image size

We compare results for the effective temperature and surface gravity of HD 4747 B to brown dwarf evolutionary models. Contours from the Saumon & Marley (2008) grids are shown in Figure 10. The range of viable T_eff values is indicated by shaded regions based on the results from spectral fitting (Table 4), where dark gray corresponds to the HK₁ best fit. The vertical extent of the plot is restricted to the range of log(g) values from Table 4. In Figure 10, this corresponds to the entire plot. Given that log(g) is only loosely constrained from spectral fitting, we also display the surface gravity uncertainty from the companion dynamical mass and model-dependent radius. As discussed in Section 4.2, the radius is estimated from the gyrochronology age using the same brown dwarf evolutionary models for self-consistency.

Zoom In Zoom Out Reset image size

Comparing contours and intersections from the Saumon & Marley (2008) grid to the allowable effective temperature and log(g) range, we find a tension (∼2σ) between the companion dynamical mass $m={65.3}_{-3.3}^{+4.4}\,{M}_{\mathrm{Jup}}$ and that predicted from evolutionary models. At an age of $\tau ={3.3}_{-1.9}^{+2.3}\,\mathrm{Gyr}$, the mass, effective temperature, and surface gravity are not coincident (Figure 9). In order to rectify this inconsistency, either (i) the dynamical mass would have to increase by ≈3 M_Jup, (ii) the effective temperature would have to decrease by ≈100 K, or (iii) the companion age would need to be lowered by several hundred megayears, which increases the companion radius and decreases log(g) by ≈0.1 dex. Each of these options can be tested with either further astrometric monitoring, higher resolution spectra obtained over a wider passband, or direct measurements of the host star radius using interferometry respectively. In the latter case, the age is inferred from the parent star based on its position in an HR diagram (Valenti & Fischer 2005; Crepp et al. 2012a). This is a model-dependent procedure that should be compared to gyrochronology in search of discrepancies (Brown 2014; Maxted et al. 2015). While age estimates are often suspect, we find that increasing the companion age only exacerbates the inconsistency, while decreasing the companion age would reduce the surface gravity and thus offer better overlap between HD 4747 B and evolutionary model contours.