# A Search for Polarized Thermal Emission from Directly Imaged Exoplanets and Brown Dwarf Companions to Nearby Stars

Rebecca Jensen-Clem, Maxwell Millar-Blanchaer, Rob van Holstein, Dimitri Mawet, James Graham, Sujan Sengupta, Mark Marley, Frans Snik, Arthur Vigan, Sasha Hinkley, et al.

<figure><img src="image_1.png" /></figure>

► To cite this version:

Rebecca Jensen-Clem, Maxwell Millar-Blanchaer, Rob van Holstein, Dimitri Mawet, James Graham, et al.. A Search for Polarized Thermal Emission from Directly Imaged Exoplanets and Brown Dwarf Companions to Nearby Stars. The Astronomical Journal, 2020, 160 (6), pp.286. 10.3847/1538-3881/abc33d. insu-03373038

HAL Id: insu-03373038

https://insu.hal.science/insu-03373038v1

Submitted on 3 Aug 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_3.png" /></figure>

# A Search for Polarized Thermal Emission from Directly Imaged Exoplanets and Brown Dwarf Companions to Nearby Stars

Rebecca Jensen-Clem¹, Maxwell A. Millar-Blanchaer², Rob G. van Holstein³,⁴, Dimitri Mawet²,⁵, James Graham⁶, Sujan Sengupta⁷, Mark S. Marley⁸, Frans Snik³, Arthur Vigan⁹, Sasha Hinkley¹⁰, Jos de Boer³, Julien H. Girard¹¹, Robert J. De Rosa⁴,¹², Brendan P. Bowler¹³, Sloane J. Wiktorowicz¹, Marshall D. Perrin¹⁴, Justin R. Crepp¹⁵, and Bruce Macintosh¹²

¹ Astronomy Department, University of California Santa Cruz, 1156 High St., Santa Cruz, CA 95064, USA; rjensenc@ucsc.edu

² Department of Astrophysics, California Institute of Technology, 1200 E. California Boulevard, Pasadena, CA 91101, USA

³ Leiden Observatory, Leiden University, P.O. Box 9513, 2300 RA Leiden, The Netherlands

⁴ European Southern Observatory, Alonso de Córdova 3107, Casilla 19001, Vitacura, Santiago, Chile

⁵ Jet Propulsion Laboratory, California Institute of Technology, 4800 Oak Grove Drive, Pasadena, CA 91109, USA

⁶ Astronomy Department, University of California, Berkeley, CA 94720, USA

⁷ Indian Institute of Astrophysics, Koramangala 2nd Block, Sarjapura Road, Bangalore 560034, India

⁸ NASA Ames Research Center, MS 245-3, Moffett Field, CA 94035, USA

⁹ Aix Marseille Univ, CNRS, CNES, LAM, Marseille, France

¹⁰ University of Exeter, Physics Department, Stocker Road, Exeter EX4 4QL, UK

¹¹ Space Telescope Science Institute, 3700 San Martin Drive, Baltimore, MD 21218, USA

¹² Kavli Institute for Particle Astrophysics and Cosmology, Stanford University, Stanford, CA 94305, USA

¹³ Department of Astronomy, The University of Texas at Austin, Austin, TX 78712, USA

¹⁴ Space Telescope Science Institute, Baltimore, MD 21218, USA

¹⁵ Department of Physics, University of Notre Dame, Notre Dame, IN 46556, USA

Received 2020 June 20; revised 2020 September 29; accepted 2020 October 19; published 2020 November 26

## Abstract

Aerosols in the atmospheres of cloudy gas giant exoplanets and brown dwarfs scatter and polarize these objects' thermal emission. If such an object has an oblate shape or nonuniform cloud distribution, the net degree of linear polarization can show an increase ranging from several tenths of a percent to a few percent. Modern high-contrast imaging polarimeters are now poised to detect such low-polarization signals, opening up a new window into the rotational velocities and cloud properties of substellar companions to nearby stars. In this paper, we present the results of a near-IR survey searching for linearly polarized thermal emission from a sample of two planetary-mass companions and five brown dwarf companions using GPI and SPHERE-IRDIS. We probe the subpercent linear polarization regime that typifies polarized free-floating brown dwarfs and place limits on each object's degree of linear polarization. We relate our upper limits on each target's degree of linear polarization to its rotation rate, and place our results in the context of rotation rates measured using high-resolution spectroscopy.

*Unified Astronomy Thesaurus concepts: Brown dwarfs (185); Polarimetry (1278); Exoplanets (498); Exoplanet atmospheres (487)*

## 1. Introduction

Observations of gas giant exoplanets and brown dwarfs indicate that clouds and hazes are ubiquitous in the atmospheres of substellar objects (e.g., Allard et al. 2001; Burgasser et al. 2002; Marley et al. 2002; Tsuji 2002; Marois et al. 2008; Barman et al. 2011). So far, this feature has largely hampered atmospheric characterization efforts: the red colors of young substellar objects are best explained by models that include thick photospheric clouds, but the radii predicted by these models are inconsistent with evolutionary models (e.g., Bowler et al. 2010; Liu et al. 2013; Tremblin et al. 2017). Such discrepancies may be due to poorly constrained cloud grain size distributions, depth variations, and horizontal structures.

However, clouds and hazes also offer observational opportunities: aerosols in the atmospheres of young, cloudy substellar objects scatter and polarize these objects' thermal emission (Sengupta & Krishan 2001). While a perfectly spherical, uniformly cloudy object has zero net polarization, an oblate shape or nonuniform cloud patterns (like Jupiter's belts and zones) can bring the net degree of linear polarization to ≲2% (de Kok et al. 2011; Marley & Sengupta 2011; Stolker et al. 2017). Since Ménard et al. (2002) reported the first

detection of linearly polarized brown dwarf radiation, several studies have sought to link polarimetric observations of substellar objects with their physical interpretations. For example, Miles-Páez et al. (2013, 2017a) find a correlation between the detection of polarized brown dwarf emission and rotation, where the rotation is based on $v \sin i$ and photometric variability measurements, respectively. These results are suggestive of oblateness-enhanced polarization. The observational case for polarization from surface inhomogeneities, however, is less developed. Observations of the photometric variability of the planetary-mass companion 2M1207 b, as well as of brown dwarfs of similar temperatures, indicate that young substellar companions' clouds are frequently patchy (Crossfield et al. 2014; Radigan 2014; Zhou et al. 2016). It has been shown by de Kok et al. (2011) that such patchy atmospheres can lead to time-varying polarimetric signals with amplitudes greater than 0.1%. Indeed, Miles-Páez et al. (2015) found observational evidence that the polarimetric variability of the cool dwarf TVLM 51346546 correlates with rotation. At the same time, detections of photometric variability and polarized emission have not been shown to be statistically correlated (e.g., Sengupta & Marley 2010; Miles-Páez et al. 2017b). Recently, Millar-Blanchaer et al. (2020) detected linearly
polarized thermal emission from the L/T transition binary brown dwarf system Luhman 16. While the polarization of Luhman 16 B could be explained either by an oblate, uniformly cloudy atmosphere or by nonuniform cloud coverage, the polarization of Luhman 16 A could only be explained by an inhomogeneous (e.g., banded) cloud structure that is constant over a single rotation period.

High-resolution spectroscopy has recently opened a new window into the rotation rates of planetary-mass companions (Snellen et al. 2014; Bryan et al. 2018; Xuan et al. 2020). Bryan et al. (2018) suggested that planets' spins are set during the late stages of accretion; subsequent theoretical work by Batygin (2018) and Ginzburg & Chiang (2020) found that these results are consistent with planetary spin regulation via magnetic torques between the planet and its circumplanetary disk. Rapidly maturing instrumentation for the high-resolution spectroscopy of faint, close-in companions (e.g., high-dispersion coronagraphy; Snellen et al. 2015; Mawet et al. 2017) will further complete this picture.

Polarimetry constrains the rotation rates of cloudy substellar objects in a way that is independent of and complementary to high-resolution spectroscopy and photometric monitoring (Sengupta & Marley 2010). For example, while both broad-band polarimetry and some approaches to high-resolution spectroscopy rely on atmospheric models to connect their observables to rotation rates, the point of connection is different. A high-resolution spectrum is cross-correlated with a model atmosphere broadened to the instrumental resolution, and the resulting cross-correlation function (CCF) is compared with model CCFs representative of a model atmosphere broadened to a given $v \sin i$ as well as the instrumental resolution (this model is also offset by a radial velocity). Alternatively, the high-resolution spectrum can be cross-correlated with spectra of slowly rotating brown dwarfs with similar spectral types to the target (Reid et al. 2002; Zapatero Osorio et al. 2006). In contrast, a measurement of an object's degree of linear polarization is related to $v \sin i$ by way of a model spectrum with a particular oblateness and set of dust scattering parameters (Sengupta & Marley 2009, 2010). Furthermore, unlike photometric or spectroscopic monitoring, a single polarimetric observation can flag a target as potentially fast-rotating if it is found to be polarized.

We also note that large discrepancies between these methods may hint at additional physical processes at play: for example, if high-resolution spectroscopic or photometric observations indicate that an exoplanet is rotating so slowly that its oblateness-induced polarization would be undetectable, but the object is nevertheless polarized, we would conclude that the polarized emission is the result of patchy clouds or a circumplanetary disk (de Kok et al. 2011; Stolker et al. 2017).

In Jensen-Clem et al. (2016), we showed that the Gemini Planet Imager (GPI; Macintosh et al. 2014) is sensitive to linear polarization fractions under 2%, using the brown dwarf companion HD 19467 B as a test case. Similarly, van Holstein et al. (2017) demonstrated that the Spectro-Polarimetric High-contrast Exoplanet REsearch instrument (SPHERE; Beuzit et al. 2019) InfraRed Dual-band Imager and Spectrograph (IRDIS; Dohlen et al. 2008) can achieve subpercent polarimetric sensitivity for substellar companions to nearby stars (these results were close to the photon noise limit). The absolute polarimetric accuracy of both instruments is estimated to be <0.1% (Wiktorowicz et al. 2014; van Holstein et al. 2020).

These capabilities led to the first direct detection of a polarized companion: CS Cha b, a ~20 $M_{Jup}$ substellar companion whose large J-band degree of linear polarization (13.7% ± 0.4%) is best explained by a disk and dust envelope (Ginski et al. 2018).

In this paper, we present the results of a polarimetric survey of two planetary-mass companions and five brown dwarf companions to nearby stars using GPI and SPHERE. Our survey is complemented by the SPHERE-IRDIS polarimetric observations of twenty substellar companions presented in R. G. T. van Holstein et al. (2020, in preparation). The structure of this paper is as follows: in Section 2, we describe our target selection process; in Section 3, we describe the observations and preliminary processing used to generate Stokes cubes from each data set; in Section 4, we retrieve the total intensity and positions of each companion; in Section 5, we measure the polarization of the parent stars; in Section 6, we determine whether polarized thermal emission has been detected and place limits on each object's linear polarization; in Section 7, we relate each object's detection limits to a model of the linear polarization fraction as a function of rotational velocity; finally, our conclusions are given in Section 8.

## 2. Target Selection

We selected our targets to cover the largest possible range of spectral types, ages, and masses of substellar companions to nearby stars whose flux would allow us to probe polarizations at the 1% level. These targets are listed in Table 1. We note that observations of brown dwarf companions offer several unique advantages compared with polarimetric observations of free-floating brown dwarfs: (1) properties such as the age and metallicity of the brown dwarf can be inferred from the bright primary star; (2) in some cases, radial velocity (RV) and/or astrometric data can be combined with orbital monitoring via imaging, to yield dynamical mass measurements of the companion; (3) the formation pathways of low-mass brown dwarf companions and massive exoplanets are potentially similar, and hence the atmospheres of brown dwarf companions are of particular interest for exoplanet science. Motivated by these advantages, we selected the age/metallicity/mass benchmark companions HR 7672 B (Liu et al. 2002; Crepp et al. 2012), HD 4747 B (Crepp et al. 2016) and HD 19467 B (Crepp et al. 2014). We further selected HR 3549 B (Mawet et al. 2015), a young brown dwarf whose spectrum has been well-characterized with SPHERE/IRDIS (Mesa et al. 2016), as well as HR 2562 B, an intermediate-age brown dwarf with hints of variability that might point toward nonuniform cloud coverage (Konopacky et al. 2016; Mesa et al. 2018). Finally, we selected two planetary-mass companions, β Pictoris b (Lagrange et al. 2010) and ROXs 42 Bb (Kraus et al. 2014), whose rotational velocities have been measured using high-resolution spectroscopy (see Section 7; Snellen et al. 2014; Bryan et al. 2018).

## 3. Observations and Data Reduction

GPI, SPHERE, and most recently SCExAO/CHARIS are the only instruments currently online that offer the combination of extreme adaptive optics and polarimetric capabilities needed to detect linearly polarized thermal emission from substellar companions to nearby stars. Observations of the targets listed in Section 2 were divided between GPI and SPHERE based on these instruments' fields of view (2"/7 × 2"/7 and 11" × 11",
<table><caption>**Table 1**
Target List</caption><thead><tr><th>Target</th><th>SpT</th><th>Age</th><th>Mass (M<sub>sun</sub>)</th><th>Distance<sup>a</sup> (pc)</th></tr></thead><tbody><tr><td>ROXs 42 Bb</td><td>M8–L0<sup>b</sup></td><td>1–5 Myr<sup>b</sup></td><td>10 ± 4 (atmospheric/evolution models)<sup>b</sup></td><td>143.59<sup>+1.54</sup><sub>-1.51</sub></td></tr><tr><td>HR 3549 B</td><td>M9–L0<sup>c</sup></td><td>100–150 Myr<sup>c</sup></td><td>40–50 (atmospheric/evolution models)<sup>c</sup></td><td>95.11<sup>+0.83</sup><sub>-0.81</sub></td></tr><tr><td>β Pic b</td><td>L2 ± 1<sup>d</sup></td><td>24 ± 3 Myr<sup>e</sup></td><td>13 ± 3 (astrometric<sup>f</sup>)</td><td>19.74<sup>+0.14</sup><sub>-0.13</sub></td></tr><tr><td>HR 7672 B</td><td>L4.5 ± 1.5<sup>g</sup></td><td>2.5 ± 1.8 Gyr<sup>h</sup></td><td>68.7<sup>+2.4</sup><sub>-3.1</sub> (dynamical<sup>h</sup>)</td><td>17.71<sup>+0.03</sup><sub>-0.02</sub></td></tr><tr><td>HD 4747 B</td><td>T1 ± 2<sup>i</sup></td><td>10.74<sup>+6.75</sup><sub>-6.87</sub> Gyr<sup>j</sup></td><td>65.3<sup>+4.4</sup><sub>-3.3</sub> (dynamical<sup>j</sup>)</td><td>18.79<sup>+0.05</sup><sub>-0.04</sub></td></tr><tr><td>HR 2562 B</td><td>T2–T3<sup>k</sup></td><td>450<sup>+300</sup><sub>-250</sub> Myr<sup>k</sup></td><td>32 ± 14 (atmospheric/evolution models)<sup>k</sup></td><td>34.01<sup>+0.05</sup><sub>-0.05</sub></td></tr><tr><td>HD 19467 B</td><td>T5–T7<sup>l</sup></td><td>10.06<sup>+1.16</sup><sub>-0.82</sub> Gyr<sup>j</sup></td><td>51.9<sup>+3.6</sup><sub>-4.3</sub> (dynamical<sup>l</sup>)</td><td>32.0<sup>+0.04</sup><sub>-0.05</sub></td></tr></tbody></table>

**Notes.**

a Bailer-Jones et al. (2018).

b Kraus et al. (2014).

c Mesa et al. (2016).

d Chilcote et al. (2017).

e Bell et al. (2015).

f Dupuy et al. (2019).

g Liu et al. (2002).

h Crepp et al. (2012).

i Crepp et al. (2018).

j Wood et al. (2019).

k Mesa et al. (2018).

l Crepp et al. (2014).

<table><caption>**Table 2**
Summary of Observations</caption><thead><tr><th>Target</th><th>Date</th><th>Instrument</th><th>Observing Mode</th><th>Exp. Time (minutes)</th><th>Parallactic Rotation (°)</th><th>Airmass</th><th>Seeing (")</th><th>Measured Stellar Q/I</th><th>Measured Stellar U/I</th></tr></thead><tbody><tr><td>β Pictoris</td><td>2015-12-24</td><td>GPI</td><td>J-pol</td><td>61.97</td><td>38.89</td><td>1.07–1.10</td><td>1.02</td><td>-0.64 ± 0.03%<sup>a</sup></td><td>0.03 ± 0.01%<sup>a</sup></td></tr><tr><td>β Pictoris</td><td>2015-12-25</td><td>GPI</td><td>J-pol</td><td>210.75</td><td>110.31</td><td>1.07–1.18</td><td>0.74</td><td>-0.24 ± 0.02%<sup>a</sup></td><td>0.03 ± 0.01%<sup>a</sup></td></tr><tr><td>ROXs 42 B</td><td>2017-03-14</td><td>SPHERE</td><td>DPI-J</td><td>139.2</td><td>8.77</td><td>1.00–1.27</td><td>0.41</td><td>1.69 ± 0.38%</td><td>0.99 ± 0.37%</td></tr><tr><td>ROXs 42 B</td><td>2017-03-15</td><td>SPHERE</td><td>DPI-J</td><td>147.73</td><td>8.92</td><td>1.00–1.38</td><td>0.47</td><td>1.74 ± 0.3%</td><td>0.94 ± 0.27%</td></tr><tr><td>HD 4747</td><td>2017-09-18</td><td>GPI</td><td>H-pol</td><td>38.83</td><td>49.1</td><td>1.01–1.02</td><td>N/A</td><td>0.01 ± 0.05%</td><td>0.02 ± 0.02%</td></tr><tr><td>HD 4747</td><td>2018-08-11</td><td>GPI</td><td>H-pol</td><td>53.9</td><td>25.17</td><td>1.01–1.07</td><td>N/A</td><td>0.04 ± 0.04%</td><td>0.08 ± 0.04%</td></tr><tr><td>HD 4747</td><td>2018-08-15</td><td>GPI</td><td>H-pol</td><td>45.63</td><td>69.78</td><td>1.01–1.02</td><td>N/A</td><td>-0.1 ± 0.03%</td><td>0.06 ± 0.02%</td></tr><tr><td>HD 19467</td><td>2017-11-26</td><td>GPI</td><td>H-pol</td><td>117.9</td><td>93.96</td><td>1.04–1.52</td><td>N/A</td><td>-0.02 ± 0.01%</td><td>0.01 ± 0.01%</td></tr><tr><td>HD 19467</td><td>2017-11-27</td><td>GPI</td><td>H-pol</td><td>31.88</td><td>21.36</td><td>1.04–1.06</td><td>N/A</td><td>-0.03 ± 0.06%</td><td>0.06 ± 0.01%</td></tr><tr><td>HR 7672</td><td>2018-05-20</td><td>GPI</td><td>H-pol</td><td>39.79</td><td>12.27</td><td>1.48–1.54</td><td>N/A</td><td>-0.05 ± 0.03%</td><td>0.01 ± 0.02%</td></tr><tr><td>HR 7672</td><td>2018-08-12</td><td>GPI</td><td>H-pol</td><td>42.36</td><td>12.78</td><td>1.48–1.55</td><td>N/A</td><td>-0.1 ± 0.02%</td><td>0.004 ± 0.004%</td></tr><tr><td>HR 2562</td><td>2019-02-21</td><td>SPHERE</td><td>DPI-H</td><td>33.33</td><td>19.37</td><td>1.23–1.25</td><td>0.93</td><td>0.08 ± 0.06%</td><td>-0.03 ± 0.06%</td></tr><tr><td>HR 2562</td><td>2019-02-25</td><td>SPHERE</td><td>DPI-H</td><td>48.27</td><td>24.08</td><td>1.26–1.37</td><td>0.91</td><td>0.07 ± 0.05%</td><td>0.06 ± 0.05%</td></tr><tr><td>HR 2562</td><td>2019-03-12</td><td>SPHERE</td><td>DPI-H</td><td>48.27</td><td>27.39</td><td>1.23–1.26</td><td>1.0</td><td>0.08 ± 0.04%</td><td>-0.01 ± 0.05%</td></tr><tr><td>HR 2562</td><td>2019-03-20</td><td>SPHERE</td><td>DPI-H</td><td>48.27</td><td>24.28</td><td>1.26–1.35</td><td>0.97</td><td>0.07 ± 0.04%</td><td>0.05 ± 0.04%</td></tr><tr><td>HR 3549</td><td>2019-03-12</td><td>SPHERE</td><td>DPI-H</td><td>35.1</td><td>18.16</td><td>1.19–1.24</td><td>0.84</td><td>0.09 ± 0.03%</td><td>0.05 ± 0.05%</td></tr><tr><td>HR 3549</td><td>2019-03-16</td><td>SPHERE</td><td>DPI-H</td><td>35.1</td><td>18.17</td><td>1.19–1.24</td><td>0.92</td><td>0.1 ± 0.05%</td><td>0.06 ± 0.04%</td></tr><tr><td>HR 3549</td><td>2019-03-23</td><td>SPHERE</td><td>DPI-H</td><td>35.1</td><td>19.63</td><td>1.18–1.20</td><td>0.96</td><td>0.02 ± 0.07%</td><td>-0.04 ± 0.05%</td></tr></tbody></table>

**Note.**

a This measurement includes the instrumental polarization as well as the stellar polarization. See Section 5 for details.

respectively) and limiting host star magnitudes for optimal adaptive optics performance ($I < 9$ and $R < 11$, respectively). The observations are summarized in Table 2. The GPI and SPHERE polarimetry modes relevant to this study are described below.

Early in our survey, we chose to observe in the *J* band, as Marley & Sengupta (2011) predict that the polarization of cloudy, oblate substellar objects will be somewhat higher in the *J* band than in the *H* band or *K* band. Later in our survey, however, it became clear that GPI’s systematics were more thoroughly characterized in the *H* band due to the large number of *H*-band polarimetric observations made as part of the GPI

Exoplanet Survey (Millar-Blanchaer et al. 2017); furthermore, both GPI and SPHERE achieve higher Strehl ratios in the *H* band. Hence, we later elected to observe in the *H* band, despite the reduction in the expected polarization signal by up to a few tenths of a percent, depending on the effective temperature and surface gravity (see Marley & Sengupta 2011 Figures 5–7).

### 3.1. GPI

In GPI’s polarimetry mode, the integral field spectrograph’s dispersing prism is replaced by a Wollaston prism, and an achromatic half-wave plate (HWP) is inserted between the
coronagraph’s focal plane and Lyot mask. A polarimetric
observing sequence consists of sets of four exposures as the
wave plate is rotated through the angles 0°, 22°, 45°, and
67°.5. GPI operates in angular differential imaging (ADI) mode,
allowing the sky to rotate with respect to the telescope pupil.
The Stokes data cube describing the astronomical polarization
is computed by inverting the Mueller matrix whose elements
represent the polarization induced by the instrument and sky
rotation. A detailed description of GPI’s polarimetry mode and
Stokes cube extraction can be found in Wiktorowicz et al.
(2014) and Perrin et al. (2015).

In this study, the Stokes data cubes were constructed from
the raw data using the publicly available GPI pipeline v1.4.0
(Maire et al. 2010; Perrin et al. 2014). Our analysis deviates
from the public pipeline’s default recipes in our treatment of
the stellar and instrumental polarization. While the public
pipeline, by default, measures and subtracts the combined
instrumental and stellar polarization from each individual
frame, we subtract only the instrumental polarization values for
each HWP position measured from the GPIES campaign data
set (Millar-Blanchaer et al. 2017). Hence, we assume that the
instrumental polarization in our data follows that of the average
GPIES campaign observation—a reasonable assumption, given
the proximity and brightness of the parent stars in our sample.
This approach allows us to measure the stellar polarization
separately from the instrumental polarization (Section 6). The
approach described above, however, cannot be applied to our
J-band β Pic b data set, as the GPIES campaign data includes
only H-band observations. Hence, we must estimate and
subtract the combined instrumental polarization and stellar
polarization for this target using the GPI pipeline.

We note that, while the crosstalk from Stokes *I* to Stokes *Q* and *U* has been characterized by Wiktorowicz et al. (2012), we do not include crosstalk in our analysis, as we do not expect that it will significantly impact the upper limits presented in Section 6.

To convert the Stokes cube data units from ADU per coadd to electrons, the cube is multiplied by the detector’s gain of 3.04e-/ ADU and the number of coadds per observation. Representative Stokes *I*, *Q*, and *U* images are shown in Figure 1.

3.2. SPHERE

We observed one planetary-mass companion and two brown
dwarf companions to nearby stars (Table 2) with the SPHERE-
IRDIS dual-beam polarimetric mode (de Boer et al. 2020; van
Holstein et al. 2020). While GPI separates a single beam into
orthogonally polarized spot pairs, SPHERE/IRDIS splits the
beam in two and uses polarizers with orthogonal transmission
axes to create two orthogonally polarized images on the
detector. An HWP is then rotated through a sequence of four
positions in order to allow the Stokes parameters to be
computed via double differencing. We construct Stokes cubes
from our raw data using the publicly available IRDIS Data
reduction for Accurate Polarimetry pipeline (IRDAP;<fnref n="16" /> van
Holstein et al. 2020). This pipeline includes a Mueller matrix
model that mitigates the effects of instrumental polarization and
crosstalk, yielding a total polarimetric accuracy of *p* ≲ 0.1%,
where *p* is the degree of linear polarization (see Section 6).
Representative Stokes *I*, *Q*, and *U* images are shown in
Figure 2.

4. Retrieval of the Total Intensity and Positions of the Companions

In order to measure the separation, position angle, and flux
of each companion in our sample, we inject negative template
PSFs into the data cubes listed in Table 2 (combining those
data sets of a given target that were acquired within one month
of each other). We use Markov Chain Monte Carlo (MCMC) to
sample the posterior distributions of the template PSFs'
separations, position angles, fluxes, and PSF shape parameters
(discussed below). Through the Vortex Image Processing (VIP;
Gomez Gonzalez et al. 2017) package interface, we use the
emcee package (Foreman-Mackey et al. 2013) affine-invariant
ensemble sampler for MCMC (Gelman & Rubin 1992).

For GPI targets, the initial template PSF is generated with the AO simulation code presented in Poyneer & Macintosh (2006). The reason for adopting this template PSF rather than making use of GPI’s satellite spots (as is done in GPI’s spectroscopy mode) is that the satellite spots are elongated in GPI’s polarimetry mode due to the full-filter bandpass. For the SPHERE targets HR 3549 B and HR 2562 B, our template PSFs are the average of the short-exposure, unocculted flux frames, obtained twice per observation sequence with a neutral density filter inserted to avoid saturation. Our final SPHERE target, ROXs 42 Bb, orbits a binary star that is barely resolved in our flux frame. To find a suitable PSF to inject into this data set, we identified 53 flux frames from the ESO archive that matched our observational setup, and chose the flux frame whose negative injection gave the best preliminary fit to the companion via VIP’s Nelder–Mead-based optimization code.

We are initially interested in measuring the position and flux
of the companions in Stokes *I* only. Hence, for each individual
exposure, we sum the two images representing the flux in the
two orthogonal directions of polarization, to create a single
total-intensity image. We then inject the negative template PSF
into the resulting total-intensity data cube. We next choose a
method for mitigating residual starlight for each data set. The
planet β Pic b is faint compared to the speckles, so we use the
VIP implementation of principal component analysis (PCA) to
subtract the speckles in a 4× full width at half maximum
(FWHM) wide annulus centered on the planet’s separation. The
remaining companions, however, are bright compared to the
speckles at the same separation. We find that ADI/PCA tends
to distort the shape of these bright objects’ PSFs, and that we
generally obtain better residuals by high-pass filtering and
median combining the intensity frames after injecting the
negative template PSF (we employ a high-pass median filter
with a box size equal to the FWHM; see Section 6). We further
find that the best-fit residuals for all objects improve if we
allow for asymmetry in the template PSF, e.g., allowing the
template PSF to stretch in the *x* and *y* directions (accomplished
by adapting the pyKLIP routine for rescaling images based on
varying wavelength), and by allowing the PSF to subsequently
rotate. We incorporate these additional PSF characteristics
unless their posterior distributions reproduce the flat prior. The
posterior distributions for flux, position, and PSF shape
parameters in addition to the best-fit residuals are shown in
Figures 17–24. Our best estimate of each object’s flux and
position is taken to be the median of the posterior distribution,
with an initial $1\sigma$ uncertainty given by the 16th and 84th
percentiles of the posterior distribution.

The analysis described above gives the GPI targets’ separations
and position angles in detector units. We use pyKLIP to translate

<footnote n="16">¹⁶ https://irdap.readthedocs.io/</footnote>
<figure><img src="image_4.png" /><figcaption>(a) β Pictoris Stokes I</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>(b) β Pictoris Stokes Q</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>(c) β Pictoris Stokes U</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>(d) HD 19467 Stokes I</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>(e) HD 19467 Stokes Q</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>(f) HD 19467 Stokes U</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>(g) HR 7672 Stokes I</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>(h) HR 7672 Stokes Q</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>(i) HR 7672 Stokes U</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>(j) HD 4747 Stokes I</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>(k) HD 4747 Stokes Q</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>(l) HD 4747 Stokes U</figcaption></figure>

**Figure 1.** Stokes cubes extracted from observations of β Pictoris (combined 2015 December 24 and 25), HD 19467 (combined 2017 November 26 and 27), HR 7672 (2018 May 20), and HD 4747 (2017 September 18). Red arrows indicate the locations of the substellar companions in the Stokes *I* frames (β Pic *b* is not visible in the raw Stokes *I* frame, but the star's debris disk is evident from the pattern in Stokes *Q* and *U*).
<figure><img src="image_16.png" /><figcaption>(g) ROXs 42 B Stokes I</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>(h) ROXs 42 B Stokes Q</figcaption></figure>

<figure><img src="image_18.png" /><figcaption>(i) ROXs 42 B Stokes U</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>(a) HR 2562 Stokes I</figcaption></figure>

<figure><img src="image_20.png" /><figcaption>(b) HR 2562 Stokes Q</figcaption></figure>

<figure><img src="image_21.png" /><figcaption>(c) HR 2562 Stokes U</figcaption></figure>

<figure><img src="image_22.png" /><figcaption>(d) HR 3549 Stokes I</figcaption></figure>

<figure><img src="image_23.png" /><figcaption>(e) HR 3549 Stokes Q</figcaption></figure>

<figure><img src="image_24.png" /><figcaption>(f) HR 3549 Stokes U</figcaption></figure>

**Figure 2.** Stokes cubes extracted from observations of ROXs 42 B (2017 March 15), HR 2562 (combined 2019 February 21, 25, 2019 March 12, and 20), and HR 3549 (combined 2019 March 12, 16, and 23). Red arrows indicate the locations of the substellar companions in the Stokes *I* frames. We note that ROXs 42 Bb is clearly detected in Stokes *Q* and *U*.

these values into units of arcseconds and degrees east of north, respectively, with the appropriate uncertainties. To the GPI objects' separation uncertainties, we add in quadrature the uncertainties due to our knowledge of GPI's platescale (14.161 ± 0.021 mas) and the location of the parent star's center (the rms of our measurements of the star's center in each exposure with the median subtracted). To the median of each object's position angle, we add GPI's offset from north (0°.36 on average), and to its uncertainty, we add in quadrature our knowledge of the direction of north (0°.12 on average) and the ratio of the uncertainty of the location of the parent star's center to the median companion-star separation. See De Rosa et al. (2015, 2020) for a discussion of the uncertainties associated with GPI's platescale and position angle—we note that these values vary in time, and hence we refer to De Rosa et al. (2020) to select the appropriate values for our observations.

To the SPHERE objects' separation uncertainties, we add in quadrature the uncertainties due to our knowledge of

SPHERE's platescale (12.255 ± 0.021 mas), a star registration error of 0.1 pixels (Vigan et al. 2015; Zurlo et al. 2016), a position angle offset error of 0°.08, and a pupil angle offset error of 0°.11 (via the SPHERE User Manual 14th public release). This method follows the approach of Nielsen et al. (2017).

To convert each object's separation to physical units, we adopt the distances and associated errors inferred from Gaia Data Release 2 by Bailer-Jones et al. (2018) and reproduced here in Table 1. A summary of our astrometric measurements is given in Table 3.

## 5. Polarimetric Analysis of the Parent Stars

Because the goal of this study is to detect polarized thermal emission from substellar companions, it is necessary to rule out sources of polarization that are external to companions' atmospheres. For example, emission from the parent star and
<table><caption>**Table 3**
Relative Astrometric Measurements</caption><thead><tr><th>Target</th><th>Instrument</th><th>Median MJD</th><th>P.A. (° on detector)</th><th>P.A. (°E. of N.)</th><th>ρ (pixels)</th><th>ρ (mas)</th><th>Proj. sep. (au)</th></tr></thead><tbody><tr><td>β Pictoris b</td><td>GPI</td><td>7381.6</td><td>214.72 ± 0.19</td><td>214.93 ± 1.38</td><td>16.48 ± 0.05</td><td>233.39 ± 5.53</td><td>4.61 ± 0.11</td></tr><tr><td>HD 4747 B</td><td>GPI</td><td>8014.7</td><td>189.9 ± 0.04</td><td>190.22 ± 0.41</td><td>41.49 ± 0.04</td><td>587.5 ± 3.97</td><td>11.04 ± 0.08</td></tr><tr><td>HD 4747 B</td><td>GPI</td><td>8341.8</td><td>193.61 ± 0.04</td><td>193.89 ± 1.21</td><td>40.77 ± 0.03</td><td>577.33 ± 12.03</td><td>10.85 ± 0.23</td></tr><tr><td>HR 7672 B</td><td>GPI</td><td>8258.8</td><td>338.58 ± 0.03</td><td>338.86 ± 0.31</td><td>30.12 ± 0.03</td><td>426.57 ± 1.97</td><td>7.56 ± 0.04</td></tr><tr><td>HR 7672 B</td><td>GPI</td><td>8342.6</td><td>336.8 ± 0.08</td><td>337.08 ± 0.26</td><td>32.43 ± 0.05</td><td>459.3 ± 1.59</td><td>8.14 ± 0.03</td></tr><tr><td>HD 19467 B</td><td>GPI</td><td>8083.7</td><td>239.02 ± 0.02</td><td>239.3 ± 0.44</td><td>115.28 ± 0.04</td><td>1632.52 ± 11.46</td><td>52.23 ± 0.37</td></tr><tr><td>ROXs 42 Bb</td><td>SPHERE</td><td>7827.8</td><td>270.25 ± 0.16</td><td>270.25 ± 0.22</td><td>95.63 ± 0.25</td><td>1171.91 ± 3.85</td><td>168.27 ± 1.87</td></tr><tr><td>HR 2562 B</td><td>SPHERE</td><td>8554.5</td><td>298.1 ± 0.25</td><td>298.1 ± 0.3</td><td>56.37 ± 0.36</td><td>690.79 ± 4.38</td><td>23.49 ± 0.15</td></tr><tr><td>HR 3549 B</td><td>SPHERE</td><td>8558.6</td><td>154.09 ± 0.11</td><td>154.09 ± 0.19</td><td>67.63 ± 0.11</td><td>828.81 ± 2.33</td><td>78.83 ± 0.71</td></tr></tbody></table>

companion alike can be polarized due to the dichroic extinction induced by interstellar dust grains that are aligned with the local magnetic field. If the companion and host star's polarizations differ, however, interstellar dust is unlikely to account for the entirety of the companion's polarization signal. It is therefore necessary to measure the polarization of each parent star.

For each parent star observed with GPI, we consider an annulus extending from just outside the focal plane mask, ~2 FWHM, to ~5 FWHM. In this region of maximum residual starlight, we can most accurately measure the stellar polarization. We recall from Section 3 that each polarimetric observing sequence consists of sets of four exposures as the HWP is rotated through the angles 0°, 22°, 45°, and 67°. From each of these exposures, the GPI pipeline constructs two images from the two orthogonally polarized channels, denoted $I^+$ and $I^-$. We compute the “normalized difference” image, given by

$$ \text{Normalized Difference} = \frac{I^+ - I^-}{I^+ + I^-}. \qquad (1) $$

We then compute the average counts inside the annulus in each normalized difference image. We use GPI's Mueller matrix model to convert these normalized differences into Stokes Q and U measurements (for a detailed description of GPI's Mueller matrix model, see Perrin et al. (2015)). These are measurements of the parent stars' polarization only, given that the instrumental polarization has already been subtracted from data used to compute the normalized difference (with the exception of the β Pic b data set, as described in Section 3). The measurement of each parent stars' polarization is given in Table 2.

To compute the uncertainties on these measurements, we calculate the standard deviation of the normalized difference values at each HWP position. These values will slightly overestimate the error on each individual normalized difference measurement at the corresponding HWP position, given that we expect the normalized difference to vary with the position angle (as any real astrophysical signal rotates with the parallactic angle). At the same time, these errors are computed from a relatively small number of frames (10 frames for most observations), which may lead to an underestimation of the true error. While imperfect, these errors represent a data-driven approach to error estimation that makes no assumptions about the source of the noise in each individual frame. We then use GPI's Mueller matrix to convert these normalized difference errors into Stokes Q and U errors; for details on this process, see Appendix B of Perrin et al. (2015). These errors are given in Table 2. When comparing our companions' polarizations to

those of their host stars (see Section 6), we assume that the stellar polarization follows a Gaussian distribution whose means are given by the Stokes Q and U values measured from the normalized differences described above and whose standard deviations follow those described in this paragraph.

For each parent star observed with SPHERE, the stellar polarization is calculated using the IRDAP package, in a procedure outlined in van Holstein et al. (2020). This is similar to our approach using the GPI data, in that the stellar polarization is computed from an annulus containing residual starlight only, and the error is calculated via the standard error on the mean of the measurements for each measurement of Stokes Q and U. The results are given in Table 2. The principal difference between the GPI and SPHERE methods for computing the stellar polarization is that the GPI method propagates the full observing sequence through the Mueller matrix, producing a single set of Stokes values and associated errors, whereas the SPHERE method computes Stokes values for each HWP cycle and finds the associated errors using the statistics of these post-Mueller matrix polarization measurements.

We note that our measurements of the stellar polarization using GPI are more precise than our measurements of the stellar polarization using SPHERE. Because GPI is located at the Cassegrain focus, its optical path is relatively static; SPHERE, however, is located at the Nasmyth focus, where the instrumental polarization and crosstalk change with time. The result is that SPHERE has somewhat larger systematic errors and uncertainties.

## 6. Polarimetric Analysis of the Substellar Companions

We now consider the polarization of the substellar companions listed in Table 1. Our methodology is based on Jensen-Clem et al. (2016), and is also adopted in R. G. T. van Holstein et al. (2020, in preparation).

For each of the substellar companions, we consider an aperture at the location of the companion and a ring of comparison apertures of the same diameter and separation from the central star. We omit those comparison apertures that fall within two apertures of the companion, to avoid contamination from the outer regions of the companion's PSF. The signal at the companion's location in the Stokes Q and U frames is the average of the difference of the aperture sum at the companion's location and the comparison aperture sums. Figure 3 illustrates how the signal in $Q/I$, $U/I$, and the degree of linear polarization ($p = \sqrt{(Q/I)^2 + (U/I)^2}$) varies with aperture size, using HD 4747 B as an example. When the PSF core of an unpolarized source contains imperfectly corrected
<figure><img src="image_25.png" /><figcaption>**Figure 3.** Q/I, U/I, and p as a function of aperture diameter for HD 4747 B. An aperture diameter of one FWHM (3.5 pixels for this GPI H-band observation) is indicated by the vertical dashed line.</figcaption></figure>

hot or cold pixels, the degree of linear polarization is artificially inflated for the smallest aperture sizes—this is apparent in Figure 3. For very large apertures, the polarization can increase due to the increasing noise contributions to the companion aperture. Across all of our targets, we find an aperture diameter of one FWHM (3.5 pixels for GPI H-band observations, 2.6 pixels for GPI J-band observations, 3.9 pixels for SPHERE H-band observations, and 3.0 pixels for SPHERE J-band observations) avoids overestimating the signal, whether due to bad pixels (for small aperture sizes) or noise (for large aperture sizes).

The signal-to-noise ratio (S/N) at the companion's location in the Stokes Q and U frames is the signal (described in the previous paragraph) divided by the standard deviation of the comparison aperture sums. As before, observations taken less than one month apart are combined to create a single Stokes cube (prior to combining data sets, we computed the S/N at the location of each companion in Stokes Q and U, and found no significant detections, with the exception of HR 7672 B, which is described in Section 6.1). We note that a measurement of the S/N at the location of β Pic b must take into account the effects of the circumstellar disk; hence, we apply a high-pass median filter with a box size equal to the FWHM of 3.5 pixels to each of that data set's polarization channels before creating the final Stokes Q and U images (we note that our flux measurement of the companion in Section 4 did not include this high-pass filtering step). With the exception of ROXs 42 Bb and HR 7672 B, we find S/N < 2 in Stokes Q and U for all observations in this study. We address ROXs 42 Bb in Section 6.2 and HR 7672 B in Section 6.1; we discuss the remaining nondetections below.

Having concluded that no polarized emission is detected from HD 4747 B, HD 19467 B, β Pic b, HR 2562 B, and HR 3549 B, we proceed to place limits on these objects' Stokes Q, U, and degree of linear polarization values. In Section 4, we obtained posterior distributions for the flux of each companion inside an FWHM-diameter aperture in Stokes I. We now seek to obtain posterior distributions of the noise at the companions' locations in Stokes Q and U. We consider again the ring of comparison apertures that was used to calculate the S/N in Stokes Q and U. Each of these comparison aperture sums represents a measurement of the background at the companion's separation. The preliminary background-corrected noise distribution at the companion's location is therefore the aperture sum at the companion's location minus each of these comparison aperture sums.

We estimate the probability density function (PDF) from this preliminary distribution using SciPy's implementation of kernel density estimation (KDE) with a Gaussian kernel and bandwidth chosen via Scott's Rule (Scott 2015). The resulting PDFs (and aperture sum histograms, for visual comparison only) are shown in Figures 4 and 5. We then draw the same number of random samples from these PDFs in Stokes Q and U as were obtained for the companion in Stokes I in Section 4 via MCMC. We note that we do not include the companion's photon noise in our Stokes Q and U PDFs—following Perrin et al. (2015) Appendix B.3, we compute the covariance of the Stokes parameters from the measured intensities of the companion and find that the companions' variances in Stokes Q and U are small compared to the variance of the comparison aperture sums (we revisit this point in Section 6.2).

The final distribution for the linear polarization fraction, p is

$$ 
\begin{aligned}
\text{PDF}(p) &= \sqrt{\left[ \frac{\text{PDF}(Q)}{\text{PDF}(I)} \right]^2 + \left[ \frac{\text{PDF}(U)}{\text{PDF}(I)} \right]^2} \\
&= \sqrt{\text{PDF}(q)^2 + \text{PDF}(u)^2},
\end{aligned} 
\quad (2) $$

where PDF(Q) and PDF(U) refer to the PDF distributions at the companions' locations in Stokes Q and U (obtained via KDE), PDF(I) refers to the posterior distribution of the companions' flux in Stokes I (obtained using negative companion injection and MCMC in Section 4), and PDF(q) and PDF(u) are the flux-normalized versions of PDF(Q) and PDF(U). The 68.27% and 99.73% confidence intervals (CIs) for q and u are obtained by stepping through the empirical cumulative distribution function (ECDF) of q or u values to find the shortest interval in q or u that encompasses a CI of 0.6827 or 0.9973 (these CIs correspond to 1σ and 3σ for a normal distribution). To place upper limits on p, a positive-definite quantity, we simply find the values of p that correspond to the ECDF value of 0.6827 or 0.9973. These limits are given for each companion in Table 4. The PDFs of p, q, and u are plotted in Figures 6 and 7. Each of our 68.27% CI upper limits (and all but one of our 99.73% CI upper limits) for p are in the subpercent regime.

## 6.1. HR 7672 B

After computing the S/N at the location of HR 7672 B as described above, we find S/N<sub>Q</sub> = 4 (and S/N<sub>U</sub> < 2.0) for both the combined 2018 data sets and the 2017 data set. Given this significant S/N at the location of the companion in Stokes Q, our task now is to decide whether the companion and stellar polarizations differ. If they are indistinguishable, we can conclude that an effect other than the intrinsic polarization of the companion is responsible for our S/N. We proceed by finding the distribution of Q and U values at the location of the companion using the approach described in Section 6. We compare these distributions with the stellar polarization calculated in Section 5 and listed in Table 2. The results are shown in Figures 8. We see that the companion and the star are indistinguishable in Stokes U, but their PDFs in Stokes Q hint at an offset between the companion and parent star. We note two important characteristics of Figure 8(a) and (c). First, in both epochs, the mean of the companion's PDF in Stokes Q is more negative than that of the host star—an encouraging sign for the interpretation that our signal is due to the intrinsic polarization of the companion. Second, in Figure 8(c), the
<figure><img src="image_26.png" /><figcaption>(a) HD 4747 B, Q (2017-09-18)</figcaption></figure>

<figure><img src="image_27.png" /><figcaption>(b) HD 4747 B, U (2017-09-18)</figcaption></figure>

<figure><img src="image_28.png" /><figcaption>(c) HD 4747 B, Q (combined 2018-08-11 and 2018-08-15)</figcaption></figure>

<figure><img src="image_29.png" /><figcaption>(d) HD 4747 B, U (combined 2018-08-11 and 2018-08-15)</figcaption></figure>

<figure><img src="image_30.png" /><figcaption>(i) HD 19467 B, Q (all observations combined)</figcaption></figure>

<figure><img src="image_31.png" /><figcaption>(j) HD 19467 B, U (all observations combined)</figcaption></figure>

<figure><img src="image_32.png" /><figcaption>(k) β Pic b, Q (all observations combined)</figcaption></figure>

<figure><img src="image_33.png" /><figcaption>(l) β Pic b, U (all observations combined)</figcaption></figure>

**Figure 4.** GPI targets: PDFs of the aperture sums at the companions' locations in Stokes Q and U, obtained via kernel density estimation, with the aperture histograms for visual comparison. HR 7672 B is omitted (see Section 6.1).

<figure><img src="image_34.png" /><figcaption>(c) HR 2562 B, Q (all observations combined)</figcaption></figure>

<figure><img src="image_35.png" /><figcaption>(d) HR 2562 B, U (all observations combined)</figcaption></figure>

<figure><img src="image_36.png" /><figcaption>(e) HR 3549 B, Q (all observations combined)</figcaption></figure>

<figure><img src="image_37.png" /><figcaption>(f) HR 3549 B, U (all observations combined)</figcaption></figure>

**Figure 5.** SPHERE targets: PDFs of the aperture sums at the companions' locations in Stokes Q and U, obtained via kernel density estimation, with the aperture histograms for visual comparison. ROXS 42 Bb is omitted (see Section 6.2).

<table><caption>**Table 4**
Nondetection Polarization Limits</caption><thead><tr><th rowspan="2">Target</th><th colspan="3">68.27% Confidence Interval</th><th colspan="3">99.73% Confidence Interval</th></tr><tr><th>p (One-sided)</th><th>q (Two-sided)</th><th>u (Two-sided)</th><th>p (One-sided)</th><th>q (Two-sided)</th><th>u (Two-sided)</th></tr></thead><tbody><tr><td>HD 4747 B</td><td>0.14%</td><td>-0.16%, -0.04%</td><td>0.00%, 0.08%</td><td>0.25%</td><td>-0.25%, 0.10%</td><td>-0.08%, 0.14%</td></tr><tr><td>HR 7672 B</td><td>0.10%</td><td>-0.10%, -0.04%</td><td>-0.00%, 0.06%</td><td>0.23%</td><td>-0.23%, 0.02%</td><td>-0.11%, 0.10%</td></tr><tr><td>HD 19467 B</td><td>0.55%</td><td>-0.16%, 0.54%</td><td>-0.35%, 0.32%</td><td>1.25%</td><td>-0.94%, 1.14%</td><td>-0.81%, 1.12%</td></tr><tr><td>β Pic b</td><td>0.18%</td><td>-0.12%, 0.12%</td><td>-0.05%, 0.16%</td><td>0.37%</td><td>-0.34%, 0.34%</td><td>-0.24%, 0.30%</td></tr><tr><td>HR 3549 B</td><td>0.06%</td><td>-0.06%, 0.01%</td><td>-0.04%, 0.04%</td><td>0.14%</td><td>-0.13%, 0.10%</td><td>-0.11%, 0.11%</td></tr><tr><td>HR 2562 B</td><td>0.26%</td><td>-0.08%, 0.23%</td><td>-0.20%, 0.10%</td><td>0.67%</td><td>-0.59%, 0.52%</td><td>-0.60%, 0.48%</td></tr><tr><td>ROXS 42 Bb</td><td>0.19%</td><td>0.11%, 0.19%</td><td>-0.10%, -0.03%</td><td>0.29%</td><td>0.01%, 0.29%</td><td>-0.17%, 0.05%</td></tr></tbody></table>

mean of the star's PDF in Stokes Q is 6.0 standard deviations from zero. However, previous polarimetric observations of HR 7672 A show no such signal: Piirola et al. (2020) measure $q = -4 \times 10^{-4}\%$ and $u = -4 \times 10^{-4}\%$ using a broad bandpass of 400–800 nm. Per Serkowski's law of interstellar

polarization, the polarization of HR 7672 A due to interstellar dust in our study's near-IR observations would be smaller still (Serkowski et al. 1975). Furthermore, the target's close proximity (17.71 pc; Bailer-Jones et al. 2018) disfavors polarization due to interstellar dust. We can therefore speculate
<figure><img src="image_38.png" /><figcaption>(a) $\beta$ Pictoris b $p$</figcaption></figure>

<figure><img src="image_39.png" /><figcaption>(b) $\beta$ Pictoris b $q$</figcaption></figure>

<figure><img src="image_40.png" /><figcaption>(c) $\beta$ Pictoris b $u$</figcaption></figure>

<figure><img src="image_41.png" /><figcaption>(d) HD 19467 B $p$</figcaption></figure>

<figure><img src="image_42.png" /><figcaption>(e) HD 19467 B $q$</figcaption></figure>

<figure><img src="image_43.png" /><figcaption>(f) HD 19467 B $u$</figcaption></figure>

<figure><img src="image_44.png" /><figcaption>(g) HR 7672 B $p$</figcaption></figure>

<figure><img src="image_45.png" /><figcaption>(h) HR 7672 B $q$</figcaption></figure>

<figure><img src="image_46.png" /><figcaption>(i) HR 7672 B $u$</figcaption></figure>

<figure><img src="image_47.png" /><figcaption>(j) HD 4747 B $p$</figcaption></figure>

<figure><img src="image_48.png" /><figcaption>(k) HD 4747 B $q$</figcaption></figure>

<figure><img src="image_49.png" /><figcaption>(l) HD 4747 B $u$</figcaption></figure>

**Figure 6.** GPI targets: PDFs of the degree of linear polarization, $Q/I$, and $U/I$ for each target with the 68.27% or 99.73% CIs indicated. Vertical dashed line indicates $Q/I = 0$ or $U/I = 0$.
<figure><img src="image_50.png" /><figcaption>(a) ROXs 42 Bb *p*</figcaption></figure>

<figure><img src="image_51.png" /><figcaption>(b) ROXs 42 Bb *q*</figcaption></figure>

<figure><img src="image_52.png" /><figcaption>(c) ROXs 42 Bb *u*</figcaption></figure>

<figure><img src="image_53.png" /><figcaption>(d) HR 3549 B *p*</figcaption></figure>

<figure><img src="image_54.png" /><figcaption>(e) HR 3549 B *q*</figcaption></figure>

<figure><img src="image_55.png" /><figcaption>(f) HR 3549 B *u*</figcaption></figure>

<figure><img src="image_56.png" /><figcaption>(g) HR 2562 B *p*</figcaption></figure>

<figure><img src="image_57.png" /><figcaption>(h) HR 2562 B *q*</figcaption></figure>

<figure><img src="image_58.png" /><figcaption>(i) HR 2562 B *u*</figcaption></figure>

**Figure 7.** SPHERE targets: PDFs of the degree of linear polarization, *Q*/*I*, and *U*/*I* for each target with the 68.27% or 99.73% CIs indicated. Vertical dashed line indicates *Q*/*I* = 0 or *U*/*I* = 0.

that our observation of the star’s nonzero polarization is due to uncorrected instrumental polarization on the order of 0.1% (Table 2).

To decide whether the companion's signal in Stokes *Q* significantly differs from that of the star (and hence differs from instrumental polarization), we mean-combine the two epochs of observation and subtract the mean-combined stellar distribution in Stokes *Q* from the mean-combined companion distribution in Stokes *Q* (we follow the same procedure for Stokes *U*). The results are shown in Figure 9. We find that the mean of the mean-combined, star-subtracted companion distribution in Stokes *Q* is located 2.1 standard deviations from zero while the mean-combined, star-subtracted companion distribution in Stokes *U* is located 0.75 standard deviations from zero.

While not sufficient to claim a detection, these results encourages follow-up observations to determine whether

HR 7672 B is intrinsically polarized. Under the assumption that we have not detected significant polarization from the companion, we can place limits on its Stokes *Q*, *U*, and * degree of linear polarization values as in Section 6. The results are plotted in Figure 6 and tabulated in Table 4.

### 6.2. ROXs 42 Bb

We now return to our analysis of ROXs 42 Bb. Using the ring of comparison apertures, as described in Section 6, we obtain S/N > 20 at the location of the companion in Stokes *Q* and *U* for both dates of observation listed in Table 2. We also find that the flux observed from the stellar host is polarized (Table 2). There are several possible explanations for the host star’s polarization: it could have a circumstellar disk; there may be magnetic star spots (e.g., Afram & Berdyugina 2015); or it
<figure><img src="image_59.png" /><figcaption>(a) The Stokes Q PDFs of HR 7672 A and B measured from observations on 2018-05-20.</figcaption></figure>

<figure><img src="image_60.png" /><figcaption>(b) The Stokes U PDFs of HR 7672 A and B measured from observations on 2018-05-20.</figcaption></figure>

<figure><img src="image_61.png" /><figcaption>(c) The Stokes Q PDFs of HR 7672 A and B measured from observations on 2018-08-12.</figcaption></figure>

<figure><img src="image_62.png" /><figcaption>(d) The Stokes U PDFs of HR 7672 A and B measured from observations on 2018-08-12.</figcaption></figure>

**Figure 8.** Comparison of the polarization of HR 7672 A and B in two epochs of observation. Vertical dashed line indicates *Q/I* = 0 or *U/I* = 0.

<figure><img src="image_63.png" /><figcaption>(a) The mean-combined, star-subtracted companion distribution in Stokes Q.</figcaption></figure>

<figure><img src="image_64.png" /><figcaption>(b) The mean-combined, star-subtracted companion distribution in Stokes U.</figcaption></figure>

**Figure 9.** Mean-combined star-subtracted HR 7672 B distributions in Stokes Q and U are consistent with zero to within 2.1 and 0.75 standard deviations, respectively.

could be polarized by the dichroic extinction induced by dust grains between the Earth and ROXs 42 B that are aligned with the local magnetic field. The spectral energy distribution of ROXs 42 B, however, shows no evidence of a circumstellar disk (Figure 10). Vrba et al. (1976) find stellar polarizations similar to that of ROXs 42 B in their observations of the Ophiuchus Cloud: Figure 11 shows their results with our observation of ROXs 42 B overplotted. While we cannot rule

out magnetic starspots, these observations of the Ophiuchus Cloud indicate that interstellar dust is a source of the star's polarization. Our task now is to decide whether the companion's polarization significantly differs from that of its host star, and hence whether we have observed polarization due to the companion itself rather than interstellar dust.

Unlike the other companions in the survey, ROXs 42 Bb is bright compared to all starlight contributions at the same
<figure><img src="image_65.png" /><figcaption>**Figure 10.** This spectral energy distribution of ROXs 42 B shows that the observed photometry is in excellent agreement with the stellar synthetic spectrum, with no evidence of a disk at infrared or (sub)mm wavelengths. Blue circles show the observed photometry. Open circles are corrected for an extinction of $A_V = 1.7$ mag (Bowler et al. 2014) using the reddening law from Fitzpatrick (1999), assuming a ratio of total to selective extinction of $R_V = 3.1$. BT-Settl atmospheric models with $T_{eff} = 3700 \pm 100$ K and $logg = 4.0$ dex are shown in gray for comparison (Baraffe et al. 2015). Choice of physical parameters is based on the spectral type of M1 ± 1 for ROXs 42 B (Bowler et al. 2014). Conversion to effective temperature is based on Herczeg & Hillenbrand (2014). Inverted triangles represent upper limits. Photometry is from UCAC4 (B and V; Zacharias et al. 2013), USNO-A (R; Monet et al. 1998), 2MASS (J, H, and Ks; Cutri et al. 2003), WISE (W1, W2, and W3; Cutri et al. 2012), Spitzer (3–24 µm; Evans et al. 2003). Upper limits at 0.88 mm and 1.3 mm are from Andrews & Williams (2007) and Wu et al. (2020), respectively.</figcaption></figure>

<figure><img src="image_66.png" /><figcaption>**Figure 11.** Reproduction of Plate VII in Vrba et al. (1976), showing stellar polarizations measured through the Ophiuchus Cloud. Length of each line is proportional to the star's degree of the linear polarization. Angle of each line indicates the degree of linear polarization (with $p = 5\%$ and $\chi = 90^\circ$ indicated in the upper left). Our result for ROXs 42 B is plotted in red.</figcaption></figure>

separation. Hence, we will use aperture photometry rather than negative companion injection to measure the flux of the companion in Stokes I. We select an aperture diameter of 5 pixels by maximizing the S/N at the location of the companion in total intensity. We find the preliminary distribution in Stokes I at the location of the planet by subtracting each comparison aperture sum from the planet's

aperture sum. We estimate the PDF from this preliminary distribution using KDE.

In calculating the PDF at the planet's location in Stokes Q and U, we must also account for the photon noise of the bright companion (PNC). Following Perrin et al. (2015), we start by considering a single measurement of Q (or U) via double differencing:

$$d_1 = I_1^+ - I_1^- = +Q + \epsilon \quad (3)$$

$$d_2 = I_2^+ - I_2^- = -Q + \epsilon, \quad (4)$$

where $I^+$ and $I^-$ are the two orthogonally polarized channels, the numerical subscripts refer to different HWP positions, and $\epsilon$ is the bias between the channels. Q is then given by

$$Q = \frac{d_1 - d_2}{2}. \quad (5)$$

The variance in Q (or U) resulting from the PNC is therefore

$$\text{var}(Q)_{\text{PNC}} = \frac{I_1^+ + I_1^- + I_2^+ + I_2^-}{4} = \text{mean}(I_{c,\text{single}}), \quad (6)$$

where $I_{c,single}$ represents the number of electrons inside the five-pixel diameter aperture at the location of the companion in one of the orthogonally polarized channels. Our final Stokes Q and U frames, however, are the average of n measurements of Q and U, where n is the number of HWP cycles. In terms of $I_c$, our measurement of the planet's signal in the final Stokes I frame, the variance at the companion's location in Q or U due to the photon noise of the companion is given by

$$\text{var}(Q) = \frac{I_c}{4n}, \quad (7)$$

where the factor of four accounts for (1) the fact that Stokes I includes twice as many frames as Q or U and (2) we are considering the signal in only one of the two orthogonally polarized channels.

To construct our final PDF in Stokes Q and U, we draw $10^6$ samples from a Gaussian distribution whose mean is the aperture sum in Stokes Q or U at the companion's location and whose variance is given by Equation (7). As in Section 6, we find the distribution of the background at the companion's separation using a ring of comparison apertures and KDE. We draw $10^6$ samples from this background distribution and subtract them from the samples drawn from the companion distribution. The resulting PDFs of Stokes Q/I and U/I at the companion's location are plotted in Figure 12. The PDFs of Stokes Q/I and U/I measured for the host star are also plotted in Figure 12 (Section 5).

We now average our PDFs of Stokes Q/I and U/I at the companion's location for the two dates of observation and subtract the similarly averaged host star PDFs. The results are plotted in Figure 14, and show that the mean-combined star-subtracted ROXs 42 Bb distributions in Stokes Q and U are separated from zero by $3.5\sigma$ and $1.7\sigma$, respectively. This difference motivates follow-up observations to determine whether the companion's polarization is indeed different from that of its host star. However, with the information at hand, we cannot conclude that this difference is significant: the companion's Stokes Q/I and U/I PDFs have different widths depending on the date of observation, and the significance of the difference between the star and companion's Stokes Q/I
<figure><img src="image_67.png" /><figcaption>(a) The Stokes Q PDFs of ROXs 42 B and Bb measured from observations on 2017-03-14.</figcaption></figure>

<figure><img src="image_68.png" /><figcaption>(b) The Stokes U PDFs of ROXs 42 B and Bb measured from observations on 2017-03-14.</figcaption></figure>

<figure><img src="image_69.png" /><figcaption>(c) The Stokes Q PDFs of ROXs 42 B and Bb measured from observations on 2017-03-15.</figcaption></figure>

<figure><img src="image_70.png" /><figcaption>(d) The Stokes U PDFs of ROXs 42 B and Bb measured from observations on 2017-03-15.</figcaption></figure>

<figure><img src="image_71.png" /><figcaption>**Figure 13.** Comparison of the degree of linear polarization (a) and angle of linear polarization (b) between the host star ROXs 42 B, stellar companion ROXs 42 Bb, and nearby background star.</figcaption></figure>

<figure><img src="image_72.png" /></figure>

PDFs is sensitive to the choice of aperture size and centroiding methods.

For the purposes of this paper, we will therefore proceed in placing an upper limit on the companion's degree of linear polarization, assuming that the host star and companion have been equally polarized by interstellar dust. It is worth noting

that the polarization of the background star located approximately 0.5 from ROXs 42 B is consistent with this interpretation: following the same analysis as described above for ROXs 42 Bb, we find that the background star is somewhat more polarized than the ROXs 42 B system (Figure 13(a); this is consistent with its greater distance, as it is likely located behind
<figure><img src="image_73.png" /><figcaption>**Figure 14.** Mean-combined star-subtracted ROXs 42 Bb distributions in Stokes Q and U are separated from zero by 3.5σ and 1.7σ, respectively. Vertical dashed line indicates $Q/I = 0$ or $U/I = 0$.</figcaption></figure>

<figure><img src="image_74.png" /><figcaption>**Figure 15.** Median oblatenesses calculated for each target using Equation (8) and the parameters in Table 1.</figcaption></figure>

the Ophiuchus Cloud) and its angle of linear polarization is similar to that of ROXs 42 B and ROXs 42 Bb (Figure 13(b)).

Continuing with the assumption that the host star and companion have been equally polarized by interstellar dust, we compute an upper limit on the intrinsic polarization of the companion by considering the host-star-subtracted distributions in Figure 14. From these distributions, we compute the 68.27% and 99.73% confidence intervals as before. These are tabulated in Table 4 and indicated on the final histograms in Figure 7.

## 7. Relating Linear Polarization to Projected Rotational Velocities

In this section, we relate the 68.27% CI upper limits on each companion's degree of linear polarization to the rotational velocity. For the purposes of this discussion, we assume an equatorial viewing angle (i.e., that $v = v \sin i$). We will estimate each object's oblateness as a function of $v$, then compute the oblateness-induced degree of linear polarization as a function of $v$.

Following Sengupta & Marley (2010) and Barnes & Fortney (2003), we find the rotation-induced oblateness by way of the

Darwin–Radau relationship:

$$f = 1 - \frac{R_p}{R_e} = \frac{\Omega^2 R}{g} \left[ \frac{5}{2} \left( 1 - \frac{3K}{2} \right)^2 + \frac{2}{5} \right]^{-1}, \quad (8)$$

where $R_p$ and $R_e$ are the companions' polar and equatorial radii, $\Omega$ is the spin angular velocity, and $g$ is the surface gravity. Finally, $K = I/(MR^2)$, where $I$ is the companion's moment of inertia, $M$ is the mass of the companion, and $R$ is the radius of the companion.

We use the Sonora Bobcat evolutionary models for brown dwarfs and self-luminous exoplanets (M. S. Marley et al. 2020, in preparation) to relate the masses and ages given in Table 1 to the radii, surface gravities, and moments of inertia needed to calculate the oblateness using Equation (8). Compared with the evolutionary models published in Saumon & Marley (2008), the Sonora Bobcat models use different atmosphere models to find the surface boundary conditions (e.g., using improved opacities of H2, CH4 and alkali resonance lines) and include metals in the interior equation of state. These differences are described in M. S. Marley et al. (2020, in preparation).

Figure 15 shows the results. We find that our companions fall into three categories: (1) the young planetary-mass objects ROXs 42 Bb and $\beta$ Pic b have the largest oblatenesses over the range of rotational velocities considered here; (2) HR 3549 B and HR 2562 B are intermediate-age, intermediate-mass brown dwarfs, and hence have intermediate oblatenesses; (3) HD 4747 B, HR 7672 B, and HD 19467B are old, higher-mass brown dwarfs, and hence have the lowest oblateness in our sample.

For each target, we then compute self-consistent, homogeneously cloudy atmospheric models for the effective temperature and log(g) values corresponding to the appropriate Sonora Bobcat evolutionary model. Previous studies of field brown dwarfs (e.g., Cushing et al. 2008; Stephens et al. 2009) have shown that typical L-type objects have cloud sedimentation efficiency $f_{sed} \sim 2$, and we adopt this value for our study. A larger value of $f_{sed}$, with thinner clouds, would produce negligible polarization, while a smaller value, with more extended, optically thick clouds, would produce somewhat greater polarization. Given the observational limits and the uncertain values for the effective temperature and log(g) of each object, this intermediate, typical, value of $f_{sed} = 2$ captures the cloud behavior without introducing an excessive number of parameters.
<figure><img src="image_75.png" /><figcaption>**Figure 16.** Degree of linear polarization corresponding to most oblate models in Figure 15 and the 68.27% CI *p* upper limits calculated for each target.</figcaption></figure>

For each of these atmospheric models, we use the multiple scattering polarization code presented in Sengupta & Marley (2009, 2010) to compute the disk-integrated linear polarization fraction for each object as a function of oblateness and rotational velocity, assuming an equatorial viewing angle (*i* = 90°)—the viewing angle that gives the highest polarization fraction for objects with an oblate shape or patchy clouds. We assume that the disk-integrated polarization is due to oblateness alone. Figure 16 plots the linear polarization fractions as a function of the rotational velocity for each object in our sample, with the exception of HR 19467 B, which, as a T-dwarf, is not compatible with our cloud models. We note that the polarization of ROXs 42 Bb is not shown for rotational velocities greater than 79 km s⁻¹, as this velocity corresponds to the planet's maximum possible oblateness value (*f* = 0.44).

Figure 16 shows our 68.27% CI *p* upper limits and linear polarization models as a function of rotational velocity. An intersection of these curves constrains the companions' projected rotational velocities, assuming *i* = 90°. Panels (b) and (c) of that figure show that our 68.27% CI *p* upper limits are insufficient to constrain the rotational velocities of the brown dwarf companions in our sample. Improving our upper limits by a factor of three to five would be sufficient to constrain these objects' rotational velocity in a future survey.

Figure 16(a) shows that our upper limit is sufficient to constrain the rotational velocity of β Pic b, but not ROXs 42 Bb. For β Pic b, we find *v* < 68 km s⁻¹. Of course, this object's rotational velocity has previously been measured using high-resolution spectroscopy: Snellen et al. (2014) measure *v* sin *i* = 25 ± 3 km s⁻¹ for β Pic b. Similarly, Bryan et al. (2018) measure *v* sin *i* = 9.5±2.1 km s⁻¹ for ROXs 42 Bb. Using the same models used to produce the linear polarization fraction versus *v* values plotted in Figure 15(a), we conclude that the measured rotational velocities of ROXs 42 Bb and β Pic b would produce *J*- and *H*-band linear polarization fractions of *p* < 10⁻³. Detecting such small polarization fractions with current instrumentation (GPI or SPHERE) would require significant improvements to our ability to subtract instrumental polarization and other instrument systematics (Wiktorowicz et al. 2014; van Holstein et al. 2020).

## 8. Conclusions and Future Work

We present a broadband near-IR polarimetric survey of two planetary-mass companions and five brown dwarf companions to nearby stars, finding that the degree of linear polarization is less

than 0.6% (68.27% CI) for all targets. We relate our upper limits to models of linear polarization as a function of rotational velocity, and place these results in the context of previous measurements of these objects' spins using high-resolution spectroscopy.

It is the combination of such diverse techniques that will allow for new avenues of exoplanet characterization in the coming years. For example, Bryan et al. (2020) combined high-resolution spectroscopy with time-resolved photometry and astrometry to constrain the true obliquity of a planetary-mass companion outside of our own Solar System for the first time. However, the main limitation of this approach was the lack of any knowledge of the sky-projected spin-orbit angle; measurements of the oblateness-induced angle of polarization would address this limitation, allowing for a much tighter constraint on the true obliquity.

As planning for the next generation of ground and space-based observatories continues, it is crucial that we design our instruments to support these varied measurements. While the James Webb Space Telescope does not have polarimetric capabilities, the Roman Space Telescope (Akeson et al. 2019) Coronagraph Instrument has an imaging polarimetry mode that, for high-S/N targets, will measure the linear polarization fraction to <3% systematic accuracy at submicron wavelengths (e.g., Bailey et al. 2018). The next generation of ground-based giant segmented mirror telescopes (GSMTs), including the Thirty Meter Telescope (TMT), Giant Magellan Telescope (GMT), and Extremely Large Telescope (ELT), do not include AO-assisted polarimetric modes with <1% accuracy in their suite of first light instruments. Second-generation GSMT instruments dedicated to high-contrast imaging (e.g., TMT/PSI; Fitzgerald et al. 2019) provide an opportunity to study a wide range of exoplanets and substellar companions in polarized light. Ongoing upgrades to the GPI, SPHERE, and SCExAO adaptive optics systems, as well as a new coronagraphic polarimetry mode planned for Keck/NIRC2 (Ragland et al. 2019), offer near-term opportunities to study faint exoplanets in polarized thermal emission in the subpercent linear polarization regime.

This material is based upon work supported by the National Science Foundation under grant No. 1816920 and 1816341. R. J.-C. acknowledges the support of the The Miller Institute for Basic Research in Science at the University of California, Berkeley. The authors thank Marta L. Bryan for many useful conversations. The research of F. Snik and J. de Boer leading to these results has received funding from the European Research Council under ERC Starting Grant agreement 678194 (FALCONER). J. R. C. acknowledges support from the NASA Early
Career Fellowship (NNX13AB03G) and NSF CAREER Fellowship (award #1654125) programs.

Facilities: Gemini(GPI), VLT(SPHERE).

Software: astropy (Astropy Collaboration et al. 2013), VIP (Gomez Gonzalez et al. 2017), GPI Data Pipeline (Perrin et al. 2014), PyKLIP (Wang et al. 2015), IRDAP (van Holstein et al. 2020), emcee (Foreman-Mackey et al. 2013) SciPy (Virtanen et al. 2020).

## Appendix

### The Results of Negative Companion Injection Coupled with MCMC

The corner plots in Figures 17–22 show the posterior distributions corresponding to the parameters of the negative companions injected into each data set. The Stokes *I* frames before and after the negative companion injection are shown in Figures 23 and 24. See Section 4 for details.

<figure><img src="image_76.png" /><figcaption>Figure 17. HD 4747 B. (a) MJD = 8014.7. (b) MJD = 8341.8.</figcaption></figure>

<figure><img src="image_77.png" /><figcaption>Figure 18. HR 7672 B. (a) MJD = 8258.8. (b) MJD = 8342.6.</figcaption></figure>
<figure><img src="image_78.png" /><figcaption>Figure 19. β Pic b.</figcaption></figure>

<figure><img src="image_79.png" /><figcaption>Figure 21. HR 3549 B.</figcaption></figure>

<figure><img src="image_80.png" /><figcaption>Figure 20. HD 19467 B.</figcaption></figure>

<figure><img src="image_81.png" /><figcaption>Figure 22. HR 2562 B.</figcaption></figure>
<figure><img src="image_82.png" /><figcaption>(a) HD 4747 B, Original, MJD = 8014.7</figcaption></figure>

<figure><img src="image_83.png" /><figcaption>(b) HD 4747 B, Subtracted, MJD = 8014.7</figcaption></figure>

<figure><img src="image_84.png" /><figcaption>(c) HD 4747 B, Original, MJD = 8341.8.</figcaption></figure>

<figure><img src="image_85.png" /><figcaption>(d) HD 4747 B, Subtracted, MJD = 8341.8.</figcaption></figure>

<figure><img src="image_86.png" /><figcaption>(e) HR 7672 B, Original, MJD = 8258.8</figcaption></figure>

<figure><img src="image_87.png" /><figcaption>(f) HR 7672 B, Subtracted, MJD = 8258.8</figcaption></figure>

<figure><img src="image_88.png" /><figcaption>(g) HR 7672 B, Original, MJD = 8342.6</figcaption></figure>

<figure><img src="image_89.png" /><figcaption>(h) HR 7672 B, Subtracted, MJD = 8342.6</figcaption></figure>

<figure><img src="image_90.png" /><figcaption>(i) HD 19467 B, Original, MJD = 8083.6</figcaption></figure>

<figure><img src="image_91.png" /><figcaption>(j) HD 19467 B, Subtracted, MJD = 8083.6</figcaption></figure>

<figure><img src="image_92.png" /><figcaption>(k) HR 3549 B, Original, MJD = 8554.6</figcaption></figure>

<figure><img src="image_93.png" /><figcaption>(l) HR 3549 B, Subtracted, MJD = 8554.6</figcaption></figure>

<figure><img src="image_94.png" /><figcaption>(m) HR 2562 B, Original, MJD = 8535.6</figcaption></figure>

<figure><img src="image_95.png" /><figcaption>(n) HR 2562 B, Subtracted, MJD = 8535.6</figcaption></figure>

**Figure 23.** High-pass filtered Stokes *I* images of the brown dwarf companions and the same images with the best-fit negative companion injected.
<figure><img src="image_96.png" /><figcaption>(c) β Pic b, Original, MJD
= 7380.7</figcaption></figure>

<figure><img src="image_97.png" /><figcaption>(d) β Pic b, Subtracted,
MJD = 7380.7</figcaption></figure>

<caption>**Figure 24.** PCA-reduced image of β Pictoris b, and the same image with best-fit negative companion injected.</caption>

## ORCID iDs

Rebecca Jensen-Clem https://orcid.org/0000-0003-0054-2953

Maxwell A. Millar-Blanchaer https://orcid.org/0000-0001-6205-9233

Rob G. van Holstein https://orcid.org/0000-0003-1520-8405

Dimitri Mawet https://orcid.org/0000-0002-8895-4735

Sujan Sengupta https://orcid.org/0000-0002-6176-3816

Mark S. Marley https://orcid.org/0000-0002-5251-2943

Frans Snik https://orcid.org/0000-0003-1946-7009

Arthur Vigan https://orcid.org/0000-0002-5902-7828

Sasha Hinkley https://orcid.org/0000-0001-8074-2562

Jos de Boer https://orcid.org/0000-0002-1028-7522

Julien H. Girard https://orcid.org/0000-0001-8627-0404

Robert J. De Rosa https://orcid.org/0000-0002-4918-0247

Brendan P. Bowler https://orcid.org/0000-0003-2649-2288

Sloane J. Wiktorowicz https://orcid.org/0000-0003-4483-5037

Marshall D. Perrin https://orcid.org/0000-0002-3191-8151

Justin R. Crepp https://orcid.org/0000-0003-0800-0593

Bruce Macintosh https://orcid.org/0000-0003-1212-7538

## References

Afram, N., & Berdyugina, S. V. 2015, A&A, 576, A34

Akeson, R., Armus, L., Bachelet, E., et al. 2019, arXiv:1902.05569

Allard, F., Hauschildt, P. H., Alexander, D. R., et al. 2001, ApJ, 556, 357

Andrews, S. M., & Williams, J. P. 2007, ApJ, 671, 1800

Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., et al. 2013, A&A, 558, A33

Bailer-Jones, C. A. L., Rybizki, J., Fouesneau, M., et al. 2018, AJ, 156, 58

Bailey, V. P., Bottom, M., Cady, E., et al. 2018, Proc. SPIE, 10698, 10698EP

Baraffe, I., Homeier, D., Allard, F., et al. 2015, A&A, 577, A42

Barman, T. S., Macintosh, B., Konopacky, Q. M., et al. 2011, ApJL, 735, L39

Barnes, J. W., & Fortney, J. J. 2003, ApJ, 588, 545

Batygin, K. 2018, AJ, 155, 178

Bell, C. P. M., Mamajek, E. E., & Naylor, T. 2015, MNRAS, 454, 593

Beuzit, J.-L., Vigan, A., Mouillet, D., et al. 2019, A&A, 631, A155

Bowler, B. P., Liu, M. C., Dupuy, T. J., et al. 2010, ApJ, 723, 850

Bowler, B. P., Liu, M. C., Kraus, A. L., et al. 2014, ApJ, 784, 65

Bryan, M. L., Benneke, B., Knutson, H. A., et al. 2018, NatAs, 2, 138

Bryan, M. L., Chiang, E., Bowler, B. P., et al. 2020, AJ, 159, 181

Burgasser, A. J., Marley, M. S., Ackerman, A. S., et al. 2002, ApJL, 571, L151

Chilcote, J., Pueyo, L., De Rosa, R. J., et al. 2017, AJ, 153, 182

Crepp, J. R., Gonzales, E. J., Bechter, E. B., et al. 2016, ApJ, 831, 136

Crepp, J. R., Johnson, J. A., Fischer, D. A., et al. 2012, ApJ, 751, 97

Crepp, J. R., Johnson, J. A., Howard, A. W., et al. 2014, ApJ, 781, 29

Crepp, J. R., Principe, D. A., Wolff, S., et al. 2018, ApJ, 853, 192

Crossfield, I. J. M., Biller, B., Schlieder, J. E., et al. 2014, Natur, 505, 654

Cushing, M. C., Marley, M. S., Saumon, D., et al. 2008, ApJ, 678, 1372

Cutri, R. M., Skrutskie, M. F., van Dyk, S., et al. 2003, 2MASS All-Sky Point Source Catalog (Washington, DC: NASA)

Cutri, R. M., Wright, E. L., Conrow, T., et al. 2012, yCat, 2311, 0

de Boer, J., Langlois, M., van Holstein, R. G., et al. 2020, A&A, 633, A63

de Kok, R. J., Stam, D. M., & Karalidi, T. 2011, ApJ, 741, 59

De Rosa, R. J., Nguyen, M. M., Chilcote, J., et al. 2020, JATIS, 6, 015006

De Rosa, R. J., Nielsen, E. L., Blunt, S. C., et al. 2015, ApJ, 814, L3

Dohlen, K., Langlois, M., Saisse, M., et al. 2008, Proc. SPIE, 7014, 70143L

Dupuy, T. J., Brandt, T. D., Kratter, K. M., et al. 2019, ApJL, 871, L4

Evans, N. J., Allen, L. E., Blake, G. A., et al. 2003, PASP, 115, 965

Fitzgerald, M., Bailey, V., Baranec, C., et al. 2019, BAAS, 51, 251

Fitzpatrick, E. L. 1999, PASP, 111, 63

Foreman-Mackey, D., Hogg, D. W., Lang, D., et al. 2013, PASP, 125, 306

Gelman, A., & Rubin, D. B. 1992, StaSc, 7, 457

Ginski, C., Benisty, M., van Holstein, R. G., et al. 2018, A&A, 616, A79

Ginzburg, S., & Chiang, E. 2020, MNRAS, 491, L34

Gomez Gonzalez, C. A., Wertz, O., Absil, O., et al. 2017, AJ, 154, 7

Herczeg, G. J., & Hillenbrand, L. A. 2014, ApJ, 786, 97

Jensen-Clem, R., Millar-Blanchaer, M., Mawet, D., et al. 2016, ApJ, 820, 111

Konopacky, Q. M., Rameau, J., Duchêne, G., et al. 2016, ApJL, 829, L4

Kraus, A. L., Ireland, M. J., Cieza, L. A., et al. 2014, ApJ, 781, 20

Lagrange, A.-M., Bonnefoy, M., Chauvin, G., et al. 2010, Sci, 329, 57

Liu, M. C., Fischer, D. A., Graham, J. R., et al. 2002, ApJ, 571, 519

Liu, M. C., Magnier, E. A., Deacon, N. R., et al. 2013, ApJL, 777, L20

Macintosh, B., Graham, J. R., Ingraham, P., et al. 2014, PNAS, 111, 12661

Maire, J., Perrin, M. D., Doyon, R., et al. 2010, Proc. SPIE, 7735, 773531

Marley, M. S., Seager, S., Saumon, D., et al. 2002, ApJ, 568, 335

Marley, M. S., & Sengupta, S. 2011, MNRAS, 417, 2874

Marois, C., Macintosh, B., Barman, T., et al. 2008, Sci, 322, 1348

Mawet, D., David, T., Bottom, M., et al. 2015, ApJ, 811, 103

Mawet, D., Ruane, G., Xuan, W., et al. 2017, ApJ, 838, 92

Ménard, F., Delfosse, X., & Monin, J.-L. 2002, A&A, 396, L35

Mesa, D., Baudino, J.-L., Charnay, B., et al. 2018, A&A, 612, A92

Mesa, D., Vigan, A., D'Orazi, V., et al. 2016, A&A, 593, A119

Miles-Páez, P. A., Pallé, E., & Zapatero Osorio, M. R. 2017a, MNRAS, 472, 2297

Miles-Páez, P. A., Zapatero Osorio, M. R., Pallé, E., et al. 2013, A&A, 556, A125

Miles-Páez, P. A., Zapatero Osorio, M. R., & Pallé, E. 2015, A&A, 580, L12

Miles-Páez, P. A., Zapatero Osorio, M. R., Pallé, E., et al. 2017b, MNRAS, 466, 3184

Millar-Blanchaer, M. A., Esposito, T. M., Stahl, K., et al. 2017, Proc. SPIE, 10407, 104070V

Millar-Blanchaer, M. A., Girard, J. H., Karalidi, T., et al. 2020, ApJ, 894, 42

Monet, D., Canzian, B., Harris, et al. 1998, yCat, 1243, 0

Nielsen, E. L., De Rosa, R. J., Rameau, J., et al. 2017, AJ, 154, 218

Perrin, M. D., Duchene, G., Millar-Blanchaer, M., et al. 2015, ApJ, 799, 182
Perrin, M. D., Maire, J., Ingraham, P., et al. 2014, Proc. SPIE, 9147, 91473J

Pirola, V., Berdyugin, A., Frisch, P. C., et al. 2020, A&A, 635, A46

Poyneer, L. A., & Macintosh, B. A. 2006, OExpr, 14, 7499

Radigan, J. 2014, ApJ, 797, 120

Ragland, S., Mawet, D., Alvarez, C., et al. 2019, in Proc. Sixth AO4ELT Conf. (Québec: AO4ELT6), http://ao4elt6.copl.ulaval.ca/proceedings.html

Reid, I. N., Kirkpatrick, J. D., Liebert, J., et al. 2002, AJ, 124, 519

Saumon, D., & Marley, M. S. 2008, ApJ, 689, 1327

Scott, D. W. 2015, Multivariate Density Estimation: Theory, Practice, and Visualization (Hoboken, NJ: Wiley)

Sengupta, S., & Krishan, V. 2001, ApJL, 561, L123

Sengupta, S., & Marley, M. S. 2009, ApJ, 707, 716

Sengupta, S., & Marley, M. S. 2010, ApJL, 722, L142

Serkowski, K., Mathewson, D. S., & Ford, V. L. 1975, ApJ, 196, 261

Snellen, I. A. G., Brandl, B. R., de Kok, R. J., et al. 2014, Natur, 509, 63

Snellen, I., de Kok, R., Birkby, J. L., et al. 2015, A&A, 576, A59

Stephens, D. C., Leggett, S. K., Cushing, M. C., et al. 2009, ApJ, 702, 154

Stolker, T., Min, M., Stam, D. M., et al. 2017, A&A, 607, A42

Tremblin, P., Chabrier, G., Baraffe, I., et al. 2017, ApJ, 850, 46

Tsuji, T. 2002, ApJ, 575, 264

van Holstein, R. G., Girard, J. H., de Boer, J., et al. 2020, A&A, 633, A64

van Holstein, R. G., Snik, F., Girard, J. H., et al. 2017, Proc. SPIE, 10400, 1040015

Vigan, A., Gry, C., Salter, G., et al. 2015, MNRAS, 454, 129

Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, NatMe, 17, 261

Vrba, F. J., Strom, S. E., & Strom, K. M. 1976, AJ, 81, 958

Wang, J. J., Ruffio, J.-B., De Rosa, R. J., et al. 2015, pyKLIP: PSF Subtraction for Exoplanets and Disks, v2.2, Astrophysics Source Code Library, ascl:1506.001

Wiktorowicz, S. J., Millar-Blanchaer, M., Perrin, M. D., et al. 2012, Proc. SPIE, 8446, 844691

Wiktorowicz, S. J., Millar-Blanchaer, M., Perrin, M. D., et al. 2014, Proc. SPIE, 9147, 914783

Wood, C. M., Boyajian, T., von Braun, K., et al. 2019, ApJ, 873, 83

Wu, Y.-L., Bowler, B. P., Sheehan, P. D., et al. 2020, AJ, 159, 229

Xuan, J. W., Bryan, M. L., Knutson, H. A., et al. 2020, AJ, 159, 97

Zacharias, N., Finch, C. T., Girard, T. M., et al. 2013, AJ, 145, 44

Zapatero Osorio, M. R., Martín, E. L., Bouy, H., et al. 2006, ApJ, 647, 1405

Zhou, Y., Apai, D., Schneider, G. H., et al. 2016, ApJ, 818, 176

Zurlo, A., Vigan, A., Galicher, R., et al. 2016, A&A, 587, A57