Cloud Atlas: Weak Color Modulations Due to Rotation in the Planetary-mass Companion GU Psc b and 11 Other Brown Dwarfs

As Figure 1 shows, there are two sources in close angular proximity of GU Psc b, leading to a second-order spectrum of a reference star and to a first-order spectrum of a galaxy (see also Figure 1 in N14). Also, there is a faint diffraction spike arises from the zeroth-order (undispersed) grism image of GU Psc. The faint spike superimposes upon part of the spectrum of GU Psc b. To mitigate this, we interpolate the flux density in the 1.17–1.19 μm region of the GU Psc b spectra (see Figure 2) to avoid possible contamination from GU Psc's zeroth-order diffraction spike.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

To assess possible contamination from the two nearby sources, we sum the measured count rates in the same row (x-axis) over three ranges of columns (i.e., the bracketed regions colored in magenta, blue, and red in Figure 1(b)), approximating the measured count rate in the J' (1.18–1.33 μm), water (1.37–1.47 μm), and H' bands (1.50–1.65 μm). The summed count rates are then plotted as a function of row number in Figure 1(a). In Figure 1(a), the summed count rate of GU Psc b spans roughly across rows 532–538 and is highlighted by the gray band. The gray region is similar to the aperture used for spectral extraction.¹⁹ We normalize the summed count rate of the three ranges of columns so that the peak count rate of GU Psc equals one.

Figure 1(a) provides an order-of-magnitude estimation of the contamination level—within the spectral aperture for GU Psc b, the contamination level is less than 10% of the GU Psc b's peak count rate in the J' and H' bands, but much higher (> 20%) in the water band. A more sophisticated analysis that fits 1D three-Moffat profiles to the summed count rate from 1.1 to 1.7 μm (column number 440–570) gives a similar estimate of the contamination level (13%; see also Appendix A). Because of the low signal-to-noise ratio and the moderate level of contamination in the water band, in this study we choose to focus on the J' band, H' band, and the integrated white (1.1–1.67 μm) light curve.

In Figure 3, the reference star's light curve shows a slight brightening trend at a subpercent level. A simple straight-line fit to the reference star's light curve gives a slope of $(9\pm 2)\times {10}^{-4}\,{\mathrm{hr}}^{-1}$, or 0.7% for an 8 hr baseline. This linear-brightening trend is possibly an HST systematics that leads to visit-long slopes (e.g., Berta et al. 2012; Wakeford et al. 2016). In the field of view no other source's light curve has a similar signal-to-noise ratio to verify the presence of this possible subpercent-level systematics. The low signal-to-noise ratio light curve of the nearby background galaxy fluctuates at about 4% from the brightest to the dimmest state. The maximum contamination level of the flux variation from the nearby sources is roughly equal to the product of the flux contamination fraction and the flux variation of nearby sources, which is about a level of 13% × 4% = 0.52%. The estimate of 0.52%, which is a generous upper limit considering the different light curve profiles of between GU Psc b and other sources, is much lower than the observed peak-to-trough flux variation of 2.7% (see Section 4.2) in the integrated white light curve of GU Psc b.

Zoom In Zoom Out Reset image size

We also check the distribution of flux density deviations from the median spectrum. We find only two data points with a 3.2σ deviation from the median spectral flux densities among 124 spectral bins in 66 exposures. Therefore, no statistically significant deviation from the median spectrum is found assuming spectral points that are independent of each other. Therefore, our careful inspection of spectral variations and reduced images confirms that the observed rotational modulation is intrinsic to GU Psc b.

We plot the HST median spectra together with the Gemini Near-InfraRed Spectrograph (GNIRS) spectrum (R ∼ 800) from N14 in the top panel of Figure 2. We also overplot the spectra of two field T dwarfs, T4.5 2MASS J05591914−1404488 (Burgasser et al. 2006a) and T3 2MASS J12095613−1004008 (Burgasser et al. 2004), which are normalized to the same J-band peak flux density as that of GU Psc b. Our HST/G141 observations provide the first flux density measurements of GU Psc b in the water-band absorption region of 1.1–1.2 μm. After the J-band flux normalization of the two field T dwarfs' spectra, GU Psc b's spectrum better matches the T4.5 spectrum at wavelengths of λ < 1.3 μm but better matches the T3 spectrum at longer wavelengths (λ > 1.3 μm).

In the bottom panel of Figure 2, we plot the binned ratio of the median spectra of the sixth and the third orbit, the orbits at which the broadband-integrated flux density reaches its maximum and minimum. The maximum/minimum flux ratio shows no strong wavelength-dependence for the rotational modulations, similar to that of T0 dwarf SIMP 0136 (Apai et al. 2013). After excluding the low signal-to-noise ratio water band (1.37–1.47 μm), the mean maximum/minimum flux ratio is about 3.0%. The fitted slope (m = $0.025\pm 0.020\,\mu {{\rm{m}}}^{-1}$) for the maximum/minimum flux ratio suggests that there is no significant wavelength dependence.

We plot the integrated white (1.1–1.67 μm) light curves of GU Psc b and those of the other two nearby sources in Figure 3. The light curve of GU Psc b manifests a sinusoidal profile with a period longer than the observation baseline. The sinusoidal pattern of GU Psc b's light curve is distinct from the almost flat light curve of the reference star and the choppy light curve of the nearby galaxy. The J'-band (1.18–1.33 μm) and H'-band (1.50–1.65 μm) light curves also show similar profiles as that of the integrated white light curve. Based on the integrated white light curve, the rotational modulation amplitude is at least 1.35%, or 2.7 ± 0.8% for the peak-to-trough flux variation (i.e., the ratio of the integrated flux median at the sixth orbit to that at the third orbit). This variability level is consistent with the previously reported marginal detection of peak-to-trough variability of 4% ± 1% at a timescale of ∼6 hr by Naud et al. (2017).

Because of the incomplete phase coverage, the fitted rotational period is degenerate with the amplitude for a sinusoidal model (see Appendix B). Therefore, we only place a lower limit of 8 hr on the rotational period, corresponding to the baseline of the HST observations.

The HST/G141 grism's wavelength coverage does not fully overlap with that of the 2MASS H band. By weighting the HST/G141 spectral variations in the HST J' (1.18–1.33 μm) and H' (1.50–1.65 μm) bands with the 2MASS J-band and H-band spectral response curves (Wright et al. 2010), we implicitly assume that modulation amplitudes ${\rm{\Delta }}J^{\prime} \,=2\mathrm{MASS}\,{\rm{\Delta }}J$ and ${\rm{\Delta }}H^{\prime} =2\mathrm{MASS}\,{\rm{\Delta }}H$.

To plot the HST/G141 spectral variations on a 2MASS CMD, we also adopt the 2MASS J- and H-band magnitudes (Cutri et al. 2003) as the mean magnitudes of the HST J'- and H'-band modulations, except for the most variable object, 2M2139.²⁰ We first scale the J-band peak of the 0.6–2.65 μm 2M2139 spectrum from the SpeX library (Burgasser 2014) to be the same as that of the averaged HST spectra. Then we use the scaled SpeX spectrum to calculate the mean magnitude in the 2MASS J and H bands during the HST observation of 2M2139.

The timescale of intrinsic color–magnitude variations from rotational modulations of brown dwarfs is typically on the order of hours. Meanwhile, the observed color–magnitude variations on several minute timescale are likely dominated by photon noise and/or systematics. To study the intrinsic color–magnitude variations across the rotational phase, for objects with long rotational periods (P > 5 hr) we bin the photometric points in the each HST orbit (total exposure time of ∼30–40 minutes) with the median value. We estimate the uncertainties of the median with the standard deviations of color and magnitude variations in each bin. These uncertainties are conservative estimates because they include the photon noise and the readout noise, variations from time-variable systematics (e.g., ramp effect), and intrinsic variability of objects. For objects with periods less than 5 hr, we choose not to bin the photometric points so that the cadence of the color–magnitude variations is less than 10% of the rotational phase.

The J − H color change due to rotational modulations in any individual object is much smaller compared to the J − H color evolution across the L/T transition. To visualize the small scale of color change and the large scale of color evolution in the CMD, we fit a straight line to the magnitude–color variation (i.e., M_J versus J − H) for the rotational modulations of each object. We then plot the fitted straight line of the magnitude–color variation on the CMD. We obtain the best-fit slope and y-intercept using an orthogonal distance regression algorithm with scipy.odr (see also Figure 4), which minimizes the orthogonal distance, weighted by both x- and y-axis uncertainties, between photometric points and the straight-line model. We use the covariance matrix of the fitted parameters to calculate the standard deviation of the fitted slope and y-intercept.

We note that this linear fit is only for illustrating the primary direction of the color changes with respect to the brightness variations. The color–magnitude variations could be nonlinear too. For example, the subpanel plots in Figure 5 suggest that a linear model may not fit well to the color–magnitude variations of GU Psc b and 2M2228. As one of the simple, periodic nonlinear models, we fit an ellipse to the color–magnitude variations and show the direction of the color change with the semimajor axis of the ellipse. The fitted color–magnitude changes of the ellipse are not unique (e.g., the elliptical and circular fit for GU Psc b in Figure 5) and are sensitive to the phase coverage of the light curves.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

For reference, we also plot the color–magnitude variations of a blackbody with varying temperatures, either as hotter or cooler from that with the effective temperature for each object as a dashed gray line in the subpanels of Figure 5. In the near-IR J and H bands, a sum of blackbodies with different temperatures is still similar to a Planck function (Schwartz & Cowan 2015). Therefore, this model acts as a toy model of a heterogeneous photosphere with a mixture of blackbodies at different temperatures.

In the left panel of Figure 5, we present the trajectories of the rotational modulations of 12 ultracool atmospheres spanning a broad spectral type range. We plot the colors and magnitudes of field brown dwarfs whose with a signal-to-noise ratio > 10 for J- and H-band photometry from the updated catalog of Dupuy & Liu (2012). The solid gray line is the empirically derived polynomial function of magnitude versus spectral type from Dupuy & Liu (2012). In the 3 × 4 subpanels in Figure 5, we show the color–magnitude variations of individual objects and the model fitting results. We plot the linear fit and the blackbody models for every object. For GU Psc b and 2M2228 that demonstrate nonlinear color–magnitude variations, we also plot the best-fit ellipses for reference. For 2M2228, the nonlinearity of variations is more apparent with the HST orbital photometric points.

The slope of the fitted lines qualitatively shows how much the J − H colors change with respect to the J-band flux variations. Among the plotted objects, 2M2139 demonstrates the largest color change (${\rm{\Delta }}(J^{\prime} -H^{\prime} )=0.0378\pm 0.0015$), but it is still much smaller than its J'-band magnitude variations (ΔJ' ∼ 0.300 ± 0.001). All of the plotted objects show less change in their colors compared to their modulation amplitudes. In other words, all objects show only weak color changes in their rotational modulations.

In the CMD, a positive slope suggests that the object becomes brighter and redder, and a negative slope suggests that the object becomes brighter and bluer. Most of the objects show negative slopes, especially among L dwarfs. Among the plotted objects, only GU Psc b and 2M2228 show positive slopes. Our fitted slope uncertainty distribution suggests that there is about a 20% probability for GU Psc b to have a negative slope as seen in other L dwarfs. Therefore, the fitted slope of GU Psc b's color–magnitude variations is not significantly different from that of the other L/T transition objects. There are two possible interpretations for the possible positive slope—either there is a phase shift between J'- and H'-band light curves of GU Psc b, as in the case of 2M2228 (phase shift of 15° ± 2°; Buenzli et al. 2012), or the H-band modulation amplitude is indeed higher than that of the J band. A longer baseline observation is needed to understand and verify if the color change of GU Psc b is indeed distinct from that of other L dwarfs.

The modulations of GU Psc b described here are the first in a planetary-mass object at the end of the L/T transition (T3-T5), confirming the previous J-band marginal detection by Naud et al. (2017). As the atmosphere evolves from L to T spectral type, the silicate cloud base presumably sank to deeper pressure (see the condensation curves in Helling & Casewell 2014; Robinson & Marley 2014), and hence the photosphere becomes less cloudy. Based on the observed rotational modulations, we argue that the photosphere of low-gravity mid-T dwarfs is not cloud free.

In addition to the modulations observed in GU Psc b, we have gathered a unique data set of rotational modulations found in planetary-mass objects across different effective temperatures, including the L7 spectral type PSOJ318 (Biller et al. 2015, 2018), T2 type SIMP0136 (Apai et al. 2013), and the T8 type Ross 458 c (Manjavacas et al. 2019). Our data set provides a useful reference point for future studies of cloud evolution across the L/T transition. Furthermore, the wavelength dependence in rotational modulations of GU Psc b will be useful in testing predictions of different cloud models on the role of gravity in shaping cloud structure, as well as on testing the hypothesis that the L/T transition, and maybe the J-band brightening too, is surface-gravity-dependent.

In any rotating atmosphere, an asymmetric brightness distribution leads to rotational modulations. The modulation amplitudes cannot be directly translated to brightness maps, as the hemisphere-integrated measurements necessarily result in information loss (although some inferences can be drawn, see Cowan & Agol 2008; Apai et al. 2013; Cowan et al. 2013; Karalidi et al. 2016). The high modulation amplitude observed in GU Psc b suggests that there are prominent rotationally asymmetric features (dark or bright) that significantly (>2%) impact even the hemisphere-integrated brightness of the photosphere. The features may be large and/or high contrast; if there are many features, these must be distributed asymmetrically. From IR rotational modulation surveys that include both high- and low-gravity brown dwarfs, we already know that brown dwarf atmospheres are ubiquitously heterogeneous (e.g., Buenzli et al. 2014; Metchev et al. 2015). Typical field brown dwarfs outside the L/T spectral type transition have peak-to-trough flux variations of less than a 2% level in the J band (Buenzli et al. 2014; Radigan et al. 2014). The observed peak-to-trough flux variations of ∼2.7% of GU Psc b is consistent with the conclusion of Radigan et al. (2014) that L/T transition objects are likely to have the largest modulation amplitudes. The high modulation amplitude of this low-gravity object provides another data point to test whether high modulation amplitude is correlated with low gravity, as claimed by Vos et al. (2017) based on a compilation of published variability amplitudes of brown dwarfs. The minimum period of 8 hr is in line with the measured timescale of periodic modulations of other brown dwarfs, ranging from as short as 1.4 hr (Buenzli et al. 2012) to 18 hr or longer (e.g., 2M2148; Metchev et al. 2015). The actual rotational modulation profile may evolve with time, as seen in long-baseline observations of multiple L/T transition brown dwarfs (e.g., Apai et al. 2017), most prominently detected in all three brown dwarfs (2M2139, 2M1324, and SIMP0136) monitored by the Spitzer Space Telescope in the Extrasolar Storms program (Yang et al. 2016; Apai et al. 2017). That study shows that the light curve evolution is the likely result of planetary-scale waves that modulated surface brightness (Apai et al. 2017), possibly through the interplay of atmospheric circulations, condensations, and cloud formation/dispersal (Tan & Showman 2017, 2019; Showman et al. 2019). These mechanisms may also be present in GU Psc b and their presence could be revealed by continuous observations over 3–4 rotational periods.

The lower limit of the rotational period of GU Psc b also provides another data point for studying the evolution of spin as a function of mass and age. Assuming that the radius of GU Psc b is roughly 1.3–1.4 Jupiter radius, as expected for a 100 Myr old and 10–12 M_Jup object (Chabrier et al. 2000), with a minimum rotational period of 8 hr, the spin velocity is about 19.6–21 kms⁻¹. A longer period (P > 8 hr) will result in a slower spin. At the age of ABDMG, which is ∼120–200 Myr, our upper limit of the spin rate of GU Psc b is similar with that of other planetary-mass companions at different ages (3–300 Myr old; e.g., Snellen et al. 2014; Zhou et al. 2016; Biller et al. 2018; Bryan et al. 2018). Because the radii contract along with the loss of interior entropy for young objects, the spin velocities increase as the objects cool with age under the conservation of angular momentum. Therefore, the spin rates of GU Psc b could reach as high as $30\,{\mathrm{kms}}^{-1}$ after radius contracting to one Jupiter radius. A better period constraint is required to test if the spin rate is consistent with the suggested universal spin–mass relation of Scholz et al. (2018; see also Zhou et al. 2019).

Both the color–magnitude variations shown in Figure 5 and the ratio of the brightest-to-dimmest spectra shown in Figure 2 demonstrate that the modulation amplitudes of GU Psc b are similar in the J' and H' bands (peak-to-trough variation of 2.6% ± 0.9% and 3.2% ± 0.9%, respectively), and hence mostly gray. Previous observations find similar modulation amplitudes in the J and H bands for L/T transition objects (e.g., Artigau et al. 2009; Radigan et al. 2012; Apai et al. 2013), although the modulation amplitudes can be different in molecular bands (e.g., water) or at longer wavelengths (e.g., K_s band; e.g., Apai et al. 2013; Biller et al. 2013). The modeling of these data suggest that spatial variations of cloud properties are responsible for the modulations (e.g., Apai et al. 2013). Detailed atmospheric modeling of space-based, high-precision time-domain data argues for simultaneous changes in cloud thickness and temperature as the cause for spatially varying cloud brightness (thin warm and thick cold clouds; Apai et al. 2013; Buenzli et al. 2015). GU Psc b also shows similarly gray rotational modulations. Therefore, we argue that correlated variations of cloud opacity and temperature are the most likely the cause of the observed gray modulations in GU Psc b.

The observed weak (non-zero) wavelength dependence may also carry information about spatial variation in the particle size distribution. A weak and positive wavelength dependence in rotational modulations could be explained by the varying presence of particles with sizes larger than one micron (e.g., see Figure 5 of Hiranaka et al. 2016). If this is true, then GU Psc b's atmosphere is in contrast to L dwarf atmospheres, which are often found to be reddened or extinguished by submicron grains (Hiranaka et al. 2016; Marocco et al. 2014; see also a retrieval analysis from Burningham et al. 2017). This would be consistent with the cloud thinning scenario that predicts larger mean particle sizes (or larger f_sed) for mid-T and smaller mean particle sizes (or smaller f_sed) for L dwarfs (Saumon & Marley 2008).

In the color–magnitude plot shown in Figure 5 even the most highly variable object (2M2139) shows only less than 20% relative difference between its J-band and the H-band modulation amplitudes (31.4 ± 0.1% versus 27.6 ± 0.1%). All other objects show similarly gray color change in modulations (${\rm{\Delta }}(J^{\prime} -H^{\prime} )\ll {\rm{\Delta }}J^{\prime}$): no known source is varying its color at a level comparable to the brightness modulations. In other words, the observed modulations in all sources are close to gray, albeit with some variety. The lack of strong color in rotational modulations in brown dwarfs with spectral types ranging from mid-L to late-T is consistent with the paradigm in which spatially heterogeneous clouds modulate the hemisphere-integrated brightness.

Given that the color changes due to atmospheric heterogeneity are different from the overall red-to-blue color evolution found across the L/T transition (a large color evolution with only a small change in absolute J-band magnitude from L8-to-T5 spectral type), we conclude that atmospheric heterogeneity alone does not directly cause the drastic color evolution across the L/T transition, at least for low-gravity atmospheres (i.e., GU Psc b, SIMP0136, PSOJ318). However, as previously suggested by Radigan et al. (2012), atmospheric heterogeneity could still affect the evolution of atmospheres over evolutionary timescales—an atmosphere with thinner clouds or more patchy cloud distribution cools more efficiently because more flux can be radiated from below the cloud base. More efficient cooling leads to a larger loss of entropy over evolutionary timescales. As a result, the loss of interior entropy is coupled with varying degrees of atmospheric heterogeneity and could still lead to the observed drastic color evolution in the L/T transition. Our observations only probe the impact of atmospheric heterogeneity with fixed interior entropy in rotational timescales. Change of cloud structure will affect the relative abundance of objects with different spectral subtypes across the L/T transition, as predicted by Saumon & Marley (2008). Observational-bias-corrected samples of brown dwarfs will be powerful to test the coupled evolution of the cloud structure and interior entropy over evolutionary timescales.

It is tempting to perceive a possible trend in Figure 5 in the slope of the magnitude–color variations from L to T dwarfs: L dwarfs become brighter and bluer, early-T dwarfs become brighter with almost no color change, and the two out of three mid-to-late T dwarfs become brighter and redder. However, we are limited by the small sample size of T dwarfs with time-resolved spectrophotometry, and the uncertain slope of GU Psc b due to incomplete phase coverage. Therefore, we remain cautious about the significance of this tentative trend. More long-term, time-resolved spectroscopy of T dwarfs is needed to verify this tentative trend and to test if there is a statistically significant difference between L and T dwarfs in the nature of their wavelength-dependent rotational modulations.