Imaging an 80 au Radius Dust Ring around the F5V
Star HD 157587

Maxwell A. Millar-Blanchaer, Jason J. Wang, Paul Kalas, James R. Graham, Gaspard Duchêne, Eric L. Nielsen, Marshall Perrin, Dae-Sik Moon, Deborah Padgett, Stanimir Metchev, et al.

► To cite this version:

Maxwell A. Millar-Blanchaer, Jason J. Wang, Paul Kalas, James R. Graham, Gaspard Duchêne, et al. Imaging an 80 au Radius Dust Ring around the F5V Star HD 157587. The Astronomical Journal, 2016, 152, 10.3847/0004-6256/152/5/128. insu-03691483

HAL Id: insu-03691483

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

Submitted on 7 Sep 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_1.png" /></figure>

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

IMAGING AN 80 au RADIUS DUST RING AROUND THE F5V STAR HD 157587

MAXWELL A. MILLAR-BLANCHAER¹,², JASON J. WANG³, PAUL KALAS³,⁴, JAMES R. GRAHAM³, GASPARD DUCHÊNE³,⁵,
ERIC L. NIELSEN⁴,⁶, MARSHALL PERRIN⁷, DAE-SIK MOON¹, DEBORAH PADGETT⁸, STANIMIR METCHEV⁹,¹⁰, S. MARK AMMONS¹¹,
VANESSA P. BAILEY⁶, TRAVIS BARMAN¹², SEBASTIAN BRUZZONE⁹, JOANNA BULGER¹³, CHRISTINE H. CHEN⁷, JEFFREY CHILCOTE²,
TARA COTTEN¹⁴, ROBERT J. DE ROSA³, RENE DOYON¹⁵, ZACHARY H. DRAPER¹⁶,¹⁷, THOMAS M. ESPOSITO³,
MICHAEL P. FITZGERALD¹⁸, KATHERINE B. FOLLETTE⁶, BENJAMIN L. GERARD¹⁶,¹⁷, ALEXANDRA Z. GREENBAUM¹⁹,
PASCALE HIBON²⁰, SASHA HINKLEY²¹, LI-WEI HUNG¹⁸, PATRICK INGRAHAM²², MARA JOHNSON-GROH¹⁶, QUINN KONOPACKY²³,
JAMES E. LARKIN¹⁸, BRUCE MACINTOSH⁶, JÉRÔME MAIRE², FRANCK MARCHIS⁴, MARK S. MARLEY²⁴, CHRISTIAN MAROIS¹⁶,¹⁷,
BRENDA C. MATTHEWS¹⁶,¹⁷, REBECCA OPPENHEIMER²⁵, DAVID PALMER¹¹, JENNIFER PATIENCE²⁶, LISA POYNEER¹¹,
LAURENT PUEYO⁷, ABHIJITH RAJAN²⁶, JULIEN RAMEAU¹⁵, FREDRIK T. RANTAKYRÖ²⁷, DMITRY SAVRANSKY²⁸,
ADAM C. SCHNEIDER²⁹, ANAND SIVARAMAKRISHNAN⁷, INSEOK SONG¹⁴, REMI SOUMMER⁷, SANDRINE THOMAS²², DAVID VEGA⁴,
J. KENT WALLACE³⁰, KIMBERLY WARD-DUONG²⁶, SLOANE WIKTOROWICZ³¹, AND SCHUYLER WOLFF¹⁹

¹Department of Astronomy and Astrophysics, University of Toronto, Toronto, ON M5S 3H4, Canada; maxmb@astro.utoronto.ca

²Dunlap Institute for Astronomy and Astrophysics, University of Toronto, Toronto, ON M5S 3H4, Canada

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

⁴SETI Institute, Carl Sagan Center, 189 Bernardo Avenue, Mountain View, CA 94043, USA

⁵Univ. Grenoble Alpes/CNRS, IPAG, F-38000 Grenoble, France

⁶Kavli Institute for Particle Astrophysics and Cosmology, Department of Physics, Stanford University, Stanford, CA, 94305, USA

⁷Space Telescope Science Institute, Baltimore, MD 21218, USA

⁸Goddard Space Flight Center, 8800 Greenbelt Road, Greenbelt, MD 20771, USA

⁹Department of Physics and Astronomy, Centre for Planetary Science and Exploration,
The University of Western Ontario, London, ON N6A 3K7, Canada

¹⁰Department of Physics and Astronomy, Stony Brook University, Stony Brook, NY 11794-3800, USA

¹¹Lawrence Livermore National Laboratory, Livermore, CA 94551, USA

¹²Lunar and Planetary Laboratory, University of Arizona, Tucson, AZ 85721, USA

¹³Subaru Telescope, NAOJ, 650 North A'ohoku Place, Hilo, HI 96720, USA

¹⁴Department of Physics and Astronomy, University of Georgia, Athens, GA 30602, USA

¹⁵Institut de Recherche sur les Exoplanètes, Département de Physique, Université de Montréal, Montréal QC H3C 3J7, Canada

¹⁶University of Victoria, 3800 Finnerty Road, Victoria, BC V8P 5C2, Canada

¹⁷National Research Council of Canada Herzberg, 5071 West Saanich Road, Victoria, BC V9E 2E7, Canada

¹⁸Department of Physics and Astronomy, University of California, Los Angeles, CA 90095, USA

¹⁹Department of Physics and Astronomy, Johns Hopkins University, Baltimore, MD 21218, USA

²⁰European Southern Observatory, Alonso de Cordova 3107, Vitacura, Santiago, Chile

²¹School of Physics, College of Engineering, Mathematics and Physical Sciences, University of Exeter, Exeter EX4 4QL, UK

²²Large Synoptic Survey Telescope, 950N Cherry Avenue, Tucson, AZ 85719, USA

²³Center for Astrophysics and Space Science, University of California San Diego, La Jolla, CA 92093, USA

²⁴Space Science Division, NASA Ames Research Center, Mail Stop 245-3, Moffett Field, CA 94035, USA

²⁵American Museum of Natural History, Department of Astrophysics, New York, NY 10024, USA

²⁶School of Earth and Space Exploration, Arizona State University, P.O. Box 871404, Tempe, AZ 85287, USA

²⁷Gemini Observatory, Casilla 603, La Serena, Chile

²⁸Sibley School of Mechanical and Aerospace Engineering, Cornell University, Ithaca, NY 14853, USA

²⁹Department of Physics and Astronomy, University of Toledo, Toledo, OH 43606, USA

³⁰Jet Propulsion Laboratory, California Institute of Technology, Pasadena, CA 91125, USA

³¹The Aerospace Corporation, 2310 E. El Segundo Boulevard, El Segundo, CA 90245, USA

Received 2016 June 9; revised 2016 August 12; accepted 2016 August 16; published 2016 October 21

ABSTRACT

We present H-band near-infrared polarimetric imaging observations of the F5V star HD 157587 obtained with the Gemini Planet Imager (GPI) that reveal the debris disk as a bright ring structure at a separation of ~80–100 au. The new GPI data complement recent *Hubble Space Telescope/STIS* observations that show the disk extending out to over 500 au. The GPI image displays a strong asymmetry along the projected minor axis as well as a fainter asymmetry along the projected major axis. We associate the minor and major axis asymmetries with polarized forward scattering and a possible stellocentric offset, respectively. To constrain the disk geometry, we fit two separate disk models to the polarized image, each using a different scattering phase function. Both models favor a disk inclination of ~70° and a 1.5 ± 0.6 au stellar offset in the plane of the sky along the projected major axis of the disk. We find that the stellar offset in the disk plane, perpendicular to the projected major axis is degenerate with the form of the scattering phase function and remains poorly constrained. The disk is not recovered in total intensity due in part to strong adaptive optics residuals, but we recover three point sources. Considering the system’s proximity to the galactic plane and the point sources’ positions relative to the disk, we consider it likely that they are background objects and unrelated to the disk’s offset from the star.

**Key words:** planet–disk interactions – stars: individual (HD 157587) – techniques: polarimetric
# 1. INTRODUCTION

Circumstellar debris disks, composed of planetesimals and dust, are remnants of the planet formation process. Therefore, their study can provide insights into the planet formation and evolution history of the systems in which they reside. The dust grain composition of a disk traces grain growth and erosion, and, if spatially resolved, disk morphology can provide evidence of dynamical interactions with nearby planets. Such an interaction can manifest as a warp (e.g., Beta Pic; Burrows et al. 1995; Mouillet et al. 1997), a stlocentric offset (e.g., HR 4796A; Wyatt et al. 1999; Telesco et al. 2000), or a sharp radial profile at the inner edge of a dust ring (e.g., Fomalhaut; Kalas et al. 2005; Quillen 2006).

Debris disks are imaged via their thermal emission in infrared or millimeter wavebands, which typically trace the location of millimeter sized bodies, or via scattered light in the visible and near-infrared (NIR), which is more sensitive to micron-sized dust. Observations of debris disks in scattered light are typically able to resolve finer spatial scales than longer wavelength observations (though ALMA's spatial resolution is now competitive), but are challenging due to the extreme contrast ratios between the faint dust-scattered light and the bright host stars. Instrumental point-spread functions (PSFs) extend the stellar emission out to angular separations where debris disks are found, obscuring the scattered light from the dust. For ground-based observations, this problem is compounded by the atmosphere, which scatters light from the PSF out to farther separations.

The Gemini Planet Imager (GPI; Macintosh et al. 2014) is an instrument on the Gemini South 8 m telescope that has been designed specifically to mitigate these challenges. It employs a high-order adaptive optics (AO) system, combined with an apodized-pupil Lyot coronagraph and an integral-field spectrograph, to image exoplanets and debris disks at angular separations down to ~0.′′1. The GPI Exoplanet Survey (GPIES) is a long-term Gemini South program targeting 600 nearby stars with the goal of discovering and characterizing young Jovian exoplanets. A secondary goal of the survey is to image and characterize debris disks. Stars with previously resolved debris disks and survey stars that exhibit infrared excesses are observed using GPI's polarimetry mode. The polarimetry mode is implemented as a rotatable half-wave plate (HWP) modulator and a Wollaston prism analyzer. This mode has been designed to take advantage of the inherent polarization of light scattered off circumstellar dust grains, to further suppress the unpolarized starlight and reveal the disk beneath (Perrin et al. 2015).

Here we present GPIES observations of the debris disk around HD 157587, an F5V star with an infrared excess $L_{IR}/L_{star} = 7.9 \times 10^{-4}$, (McDonald et al. 2012) at a distance of 107.4 pc (van Leeuwen 2007). *Hubble Space Telescope (HST)/STIS coronagraphic imaging* (GO-12998; PI Padgett) first revealed the dust-scattered light extending to a $>7′′$ radius, with a morphology resembling a fan (such as for HD 15745; Kalas et al. 2007), where the straight edge of the fan lies along the southwestern side of the nebulosity (Padgett & Stapelfeldt 2015). The inner working angle of these data corresponds to a projected separation of ~100 au. Our new scattered-light images, obtained as part of the GPIES campaign, detect the structure of the circumstellar dust in the projected 30–130 au radial region.

# 2. OBSERVATIONS AND DATA REDUCTION

We observed HD 157587 with GPI's polarimetry mode in the H-band on 2015 August 28 UT. The observations consisted of twenty-eight 90 s frames, with the HWP position angle cycling between 0°, 22°, 45°, and 67°. Throughout the sequence, the field rotated by a total of 46°. The average airmass was 1.02 and the seeing as measured by the Gemini Differential Image Motion Monitor and Multi-Aperture Scintillation Sensor was 0.′′61 and 0.′′63, respectively. The AO system telemetry reported a post-correction wavefront rms error of 216 ± 20 nm across the sequence.

The data were reduced using the GPI data reduction pipeline version 1.3 (Maire et al. 2012; Perrin et al. 2014). The raw data were dark subtracted, cleaned of correlated detector noise, bad pixel corrected, flexure corrected, and then combined into a polarization datacube (where the third dimension holds two orthogonal polarization states). Each datacube was divided by a polarized flat field and corrected for non-common path errors via a double differencing algorithm (Perrin et al. 2015). The instrumental polarization was determined by estimating the apparent stellar polarization in each polarization datacube by measuring the mean normalized difference of pixels with separations between 7 and 13 pixels from the star's location (determined from the satellite spots using a radon-transform-based algorithm; Wang et al. 2014). The estimated instrumental polarization was then subtracted from each pixel, scaled by the pixel's total intensity (Millar-Blanchaer et al. 2015). The region selected to measure the instrumental polarization was just outside of the coronagraph edge where the residual PSF flux, and hence the flux from instrumental polarization, is maximized. We assume that this area is devoid of any significant polarized structure and that any measured difference between the two polarization states is due to the instrumental polarization.

The datacubes were corrected for geometric distortion, smoothed with a Gaussian kernel ($\sigma = 1$ pixel) and then combined into a Stokes datacube by solving a set of equations that describe the linear polarization states measured in each of the individual exposures given the waveplate and sky rotation angles (Perrin et al. 2015). The Stokes datacube was subsequently converted to the radial Stokes convention ($[I, Q, U, V] \rightarrow [I, Q_r, U_r, V]$; Schmid et al. 2006). The sign convention is such that a positive $Q_r$ corresponds to a polarized intensity whose vectors are oriented perpendicular to a line connecting a given pixel to the central star and negative values are parallel to the line. Under this convention (and the assumption of low optical depth) the $U_r$ image should contain no disk flux and will only contain noise. Thus the $Q_r$ image should contain all of the disk polarized intensity as positive values. Finally, the flux of the four satellite spots was measured and flux calibration was carried out as described in Hung et al. (2015). The final $Q_r$ and $U_r$ images can be seen in Figure 1.

The polarization datacubes were also processed separately using the *pyKLIP* (Wang et al. 2015) implementation of the Karhunen–Loève Image Projection (KLIP) algorithm (Soummer et al. 2012) to attempt to recover the disk in total intensity and search for point sources (Section 3.2).
<figure><img src="image_3.png" /></figure>

**Figure 1.** Left: GPI H-band radial polarized intensity image ($Q_r$) of the HD 157587 debris disk. The image has been smoothed with a Gaussian kernel ($\sigma$ = 1 pixel). The blue circle indicates the size of the central focal plane mask (0'' 12 radius). The red circles denote the point sources seen in total intensity in Figure 4 (Section 3.2). Right: the $U_r$ image of HD 157587. For single scattering from circumstellar material, we expect no contributions to the $U_r$ image. Thus, this image can be used as a noise map for the $Q_r$ image. The image appears to be largely free of correlated structure, except at small inner working angles. Both images have been cropped from the full GPI field of view to display only the inner 2'' 6 $\times$ 2'' 6 region. No polarized emission was seen outside the cropped region.

## 3. RESULTS

### 3.1. Polarized Intensity Image

The $Q_r$ image displays an inclined, ring-like structure with a strong brightness asymmetry in the NE–SW direction (the projected minor axis). The inner edge of the ring has projected semimajor and semiminor axes of ~0'' 65 and ~0'' 2, respectively. The region interior to the inner edge of the ring appears to be cleared of any scattering material. However, the residual systematics at smaller separations in the $Q_r$ and $U_r$ images do not exclude additional dust at these smaller radii. Outside of the ring, the surface brightness decreases quickly and reaches the noise floor within our field of view (FOV), which extends to a radius of 1'' 7 along the semimajor axis of the ring.

A comparison between the $Q_r$ and $U_r$ images indicates that the ring detection is robust and that the morphology is not due to instrumental effects. The $U_r$ image appears to be dominated by uncorrected systematics interior to ~0'' 275, a region that intersects the ring near its minor axis. Outside of ~0'' 275, there appear to be no coherent structures in the $U_r$ image.

The strong NE–SW asymmetry seen in the polarized images is reminiscent of the asymmetries seen in other disks recently imaged in polarized light by GPI, for example HR 4796A, Perrin et al. (2015); HD 106906, Kalas et al. (2015); HD 131835, Hung et al. (2015); and HD 61005 (Esposito et al., 2016). In all of these disks, this asymmetry is interpreted as the disk being tilted such that the brighter side is closer to the observer and the observed brightness asymmetry is mostly due to strong forward scattering in the polarized scattering phase function. Indeed, a recent analysis of Cassini observations (albeit total intensity visible light observations) of Saturn's G and D rings indicates that collisionally generated dust is expected to be strongly forward scattering (Hedman & Stark 2015).

In addition to the NE–SW asymmetry, the $Q_r$ image also displays a mild brightness asymmetry between the SE and NW sides of the disk, visible as two main features: (a) the SE ansa appears brighter and reveals more of the backside of the disk

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

**Figure 2.** Signal-to-noise map of the GPI H-band $Q_r$ image. The ansae are clearly detected here with an S/N greater than five, while the region near the minor axis is closer to an S/N of three. In this map, the projected minor axis (in the NE direction) appears at a low S/N even though it is the brightest region in the $Q_r$ image. This is due to the elevated noise in the inner regions of the $U_r$ image that overlap with the edge of the disk.

than the NW ansa, causing the ansa to appear hook-like, and (b) the SE side of the disk appears brighter along the bright NE edge of the disk, about the NE semiminor axis. These features are confirmed in a signal-to-noise ratio (S/N) map (Figure 2), created by dividing the $Q_r$ image at each point by the standard deviation of an annulus in the $U_r$ image at the same angular separation.

To estimate the magnitude of this asymmetry, we created a custom-shaped aperture for each ansa (Figure 3). The two apertures are mirror images of each other, with the axis of symmetry coincident with the projected minor axis (a position angle of 37°, Section 4). By summing the flux in each aperture, we find the brightness ratio between the SE ansa and the NW
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Regions used to calculate the magnitude of the brightness asymmetry between the two ansae overplotted on top of the *Q*<sub>r</sub> image from Figure 1.</figcaption></figure>

ansa to be $1.15 \pm 0.02$, where the errors at each pixel are calculated in the same manner as when creating the S/N map.

This brightness asymmetry may be explained by a stellar offset, which in turn may be caused by a perturbing planet in an eccentric orbit that imparts a forced eccentricity to the dust's parent bodies (e.g., Wyatt et al. 1999). For small eccentricities, the morphology of the disk remains axisymmetric (to first order), but the host star is no longer located at the geometric center of the disk. A brightness asymmetry can then be seen as a result of one side of the disk being closer to the star than the other and receiving increased stellar irradiation. This offset also warms the closer dust gains, an effect known as pericenter glow, which can be observed in thermal emission (Wyatt et al. 1999).

## 3.2. Total Intensity Image

Each individual polarization datacube was summed across its two polarization channels to create a total intensity datacube. The entire set was then processed with *pyKLIP* using a large range of Karhunen–Loève (KL) modes, exclusion criteria and optimization regions. No disk emission was detected in any reduction. This is due in part to prominent stellar residuals resulting from imperfect AO correction caused by strong winds throughout the observation set. These winds held a roughly constant position angle during the observations, which caused the AO residuals to rotate relative to the instrument frame, mimicking the rotation of an astrophysical source. Aggressive PSF subtraction is able to suppress nearly all of this signal, but also suppresses the disk emission. In addition to the wind residuals, any angular differential imaging-based method will be subject to self-subtraction for such an azimuthally extended disk, compounding the difficulties in detecting the disk in total intensity.

Although no dust-scattered light was detected, we recovered three possible point sources in the PSF-subtracted data (Figure 4). Their measured properties are summarized in Table 1, where we have labeled them candidate companions

<figure><img src="image_6.png" /><figcaption>**Figure 4.** PSF-subtracted total intensity (Stokes *I*) image of HD 157587 at *H* band. No disk is detected, in part because of strong AO residuals caused by winds. Residuals due to the winds can be seen as nearly vertical dark streaks to the north and south of the obscured star. The three point-sources described in Table 1 are marked by red circles. We consider the region inside the green dashed lines to be dominated by wind residuals.</figcaption></figure>

(cc) 1, 2, and 3. The brightest of the sources, cc1, is confidently detected. Due to its position on the edge of the field of view, cc2 is just below a $5\sigma$ detection. The faintest of the three sources (source cc3) is recovered at less than $3\sigma$ significance. However, we find that both cc2 and cc3 are stable as a function of KL modes and appear as point sources in both our most conservative (i.e., with a low number of KL modes) and aggressive (i.e., with a high number of KL modes) reductions, which does not hold true for other low significance point-source candidates in the data. Thus, we have decided to report both alongside the one confident detection.

The flux and position of the point sources were calculated using a Gaussian matched filter. The S/N was determined by comparing the flux of the point sources with the noise at the same radial separation. Because the point sources lie outside of the region with strong wind residuals, we used a parallelogram-shaped region to mask out the wind residuals when estimating the noise (see Figure 4). To correct for algorithm throughput and to characterize the uncertainties, artificial point sources of known brightness and position were injected into the data at similar separations but at different azimuthal positions with respect to the point sources cc1, cc2, and cc3, avoiding the region with strong wind residuals. Algorithm throughput was estimated by measuring the flux of the artificial point sources after PSF subtraction. The scatter in the position and flux of the artificial planets were used as the uncertainties on the position and flux of the point sources, respectively. To obtain the total error in the astrometry, we used the reported plate scale and north angle from De Rosa et al. (2015) and added the uncertainties in quadrature. For our flux conversion, we used the flux of the satellite spots to convert the flux of the point sources to contrast units, using the standard GPI calibrations for the flux ratio of the satellite stars relative to the central PSF
<table><caption>**Table 1**
Candidate Point Source Properties</caption><thead><tr><th>Label</th><th>S/N</th><th>Separation</th><th>Position Angle</th><th>H-band Flux Ratio<sup>a</sup></th><th>Radial Separation<sup>b</sup></th></tr></thead><tbody><tr><td>cc1</td><td>6.6</td><td>1"180 ± 0"002</td><td>228°9 ± 0°2</td><td>(3.2 ± 0.8) × 10<sup>-6</sup></td><td>364 au</td></tr><tr><td>cc2</td><td>4.2</td><td>1"248 ± 0"005</td><td>195°9 ± 0°2</td><td>(2.7 ± 0.7) × 10<sup>-6</sup></td><td>380 au</td></tr><tr><td>cc3</td><td>2.8</td><td>1"002 ± 0"004</td><td>269°3 ± 0°3</td><td>(1.9 ± 0.7) × 10<sup>-6</sup></td><td>210 au</td></tr></tbody></table>

**Notes.**

<sup>a</sup> Between the source and the star.

<sup>b</sup> If on a circular orbit in the disk plane.

<figure><img src="image_7.png" /><figcaption>**Figure 5.** Sensitivity to point sources in the total intensity image. The dashed and solid blue lines are the sensitivity curves inside and outside of the wind residuals, respectively. The candidate point sources are plotted in orange. Consistent with their reported S/Ns in Table 1, only b is detected above a significance of 5σ.</figcaption></figure>

(Wang et al. 2014). The scatter in the satellite spot fluxes was used as the uncertainty in the flux conversion factor.

While it is possible that one or more of these three point sources is associated with HD 157587, we note that the star's projected position on the sky is near the galactic plane ([l, b] = [6°0, 9°4]) and it is likely that most, if not all, of these sources are background objects. We further discuss the potential relationship of these point sources to the debris disk in Section 5.

To assess our sensitivity to additional point sources in the image, we generated contrast curves for our total intensity image (Figure 5). Because of the strong wind residuals, we computed contrast curves both inside and outside of the wind residuals using the aforementioned parallelogram shape (seen in Figure 4) to define the two regions. Injected point sources were again used to correct for algorithm throughput and the reported 5σ contrasts were corrected for small number statistics (Mawet et al. 2014).

## 4. DISK MODELING

To recover basic geometric properties of the disk, we modeled the $Q_r$ image using two modified versions of the disk model presented in Millar-Blanchaer et al. (2015). The original model describes the three-dimensional dust density as a radial power law centered on the host star with a Gaussian height profile and constant aspect ratio (the ratio of the disk scale height to the radial separation). Optically thin (single) scattering is assumed and a Henyey-Greenstein (HG) function is used as a polarized scattering phase function. A disk image is

calculated by combining the dust density profile with the scattering phase function and integrating along the line of sight. The model includes nine free parameters: inner radius, $R_1$; outer radius, $R_2$; power-law index for the radial dust distribution, $\beta$; disk aspect ratio, $h_0$; average scattering cosine, $g$; inclination, $\phi$; position angle on the sky, $\theta_{PA}$; a flux normalization factor $N_0$ and a constant offset term applied to the entire image $I_0$.

Motivated by the SE-NW asymmetries seen in Figure 1, we adapted the disk model to allow for a stellar offset in the plane of the disk. This was implemented as two new free parameters: $\Delta X_1$, an offset along the projected semimajor axis of the disk (i.e., parallel to the position angle in the plane of the sky; positive offset is in the SE direction), and $\Delta X_2$, an offset in the disk plane, perpendicular to $\Delta X_1$ (positive offset is toward the backside of the disk). The true orientation of $\Delta X_2$ relative to the observer is parallel to the inclination vector $\phi$ and will have components both along the projected semiminor axis of the disk and along the observer's line of sight. For example, if the disk is edge-on ($\phi = 90°$), then $\Delta X_2$ is parallel to the line of sight and if the disk is face-on, then $\Delta X_2$ will be in the plane of the sky.

We considered two different options for the polarized scattering phase function. We first used the HG function, as in the original model, and we also considered a model where the polarized scattering phase function is described by a Rayleigh scattering function multiplied by an HG function. This second model was considered because the Rayleigh scattering function produces a peak in the polarization fraction at 90° scattering angles and therefore places increased importance on the ansae, where the hook-feature is seen. In this second model the HG function describes the scattering phase function of the (undetected) total intensity.

A similar strategy was used by Graham et al. (2007), who simultaneously fit total intensity and polarization fraction observations of the AU Mic debris disk with the same combination of an HG function and Rayleigh scattering function. One difference between their analysis and ours is that their fit includes the maximum polarization fraction, $p_{max}$, as a free parameter. In our case, with no total intensity detection, we fit only the polarized intensity, not the polarization fraction, and the maximum polarization fraction is folded into the flux normalization term, $N_0$.

Following Millar-Blanchaer et al. (2015), we fit the two disk models to the image using the affine invariant Markov Chain Monte Carlo (MCMC) ensemble sampling package `emcee` (Foreman-Mackey et al. 2013). The package uses an ensemble of “walkers,” each with its own MCMC chain, where each step of each walker depends on the state of the other walkers. Before fitting, we first applied 3 × 3 pixel binning to both the $Q_r$ image and the $U_r$ image. At H-band, the diffraction limit of
<figure><img src="image_8.png" /><figcaption>**Figure 6.** Left: images of the best-fit *H*-band, *Q*<sub>r</sub> models described in Table 2. Right: the residuals of each model subtracted off the *Q*<sub>r</sub> image. The top row displays the disk image and residuals for the model that uses only an HG function to describe the polarized scattering phase function. The bottom row displays the disk image and residuals for the model that uses an HG function combined with Rayleigh scattering. The images are displayed at the same color scale and orientation as Figure 1.</figcaption></figure>

Gemini is about 3 GPI pixels. Therefore, by binning the data, we can improve S/N without any significant loss of spatial information. This has the added benefit that each binned pixel can be considered nearly independent, thus correlated errors between the binned pixels are significantly reduced relative to the full-resolution pixels. The error estimates for each pixel were calculated as in Section 2, by taking the standard deviation of concentric annuli in the (binned) *U*<sub>r</sub> image. Preliminary fitting runs revealed that there was no evidence for an offset term, *I*<sub>0</sub> and as a result we opted to drop the offset as a free parameter.

For the first model (HG-only), the fitting procedure was run with 240 walkers and a burn-in stage of 500 iterations followed by a full run of 1300 iterations. The burn-in serves to initialize the state of the walkers and only walker positions from the full run are used to calculate the posterior. The fitting code was executed on the Edison supercomputer at the National Energy Research Scientific Computing Center (NERSC), employing one process for each walker. For each iteration, a new model was generated at the same resolution as a full GPI frame, smoothed with a Gaussian filter (σ = 1 pixel) and then binned 3 × 3, to replicate the final steps of the data reduction.

After the run, the maximum autocorrelation across all parameters was found to be 65 iterations, indicating that the chains had iterated longer than the required *O*(10 × *t*<sub>autocorrelation</sub>) for convergence (see Foreman-Mackey et al. 2013). As a second check, the ensemble chains were examined by eye after the burn-in stage and appeared to have converged. The entire procedure took 46 hr and 30 minutes to complete using 240 Intel Ivy-Bridge 2.4 GHz cores (10 nodes each with 24 cores).

We fit the data to the second model (HG + Rayleigh), using the same *emcee* setup, but only 1000 iterations were completed due to time constraints on the NERSC Edison supercomputer. For this run, the autocorrelation time was found to be 91, again, indicating that the chains should have converged.

Marginalized posterior distributions for each parameter of each of the two models were obtained by sampling the MCMC chains of each walker at intervals equal to one autocorrelation time. A summary of the results of the fitting procedure can be found in Table 2, which reports the best-fit value for each parameter as the median of its marginalized posterior, and the errors are taken to be the 68% confidence intervals. The normalization factor *N*<sub>0</sub> is considered a nuisance parameter and
<table><caption>**Table 2**
Disk Model Parameters</caption><thead><tr><th rowspan="2">Parameter</th><th rowspan="2">Symbol</th><th rowspan="2">Units</th><th rowspan="2">Prior Range</th><th colspan="2">Best-fit Values</th></tr><tr><th>HG Model</th><th>HG + Rayleigh Model</th></tr></thead><tbody><tr><td>Inner radius</td><td>R<sub>1</sub></td><td>au</td><td>40–120</td><td>81 ± 2</td><td>78.9 ± 0.8</td></tr><tr><td>Outer radius</td><td>R<sub>2</sub></td><td>au</td><td>R<sub>1</sub>–1000</td><td>216 ± 16</td><td>211<sup>+21</sup><sub>-15</sub></td></tr><tr><td>Density power-law index</td><td>β</td><td>...</td><td>0.1–4</td><td>2.23 ± 0.15</td><td>2.2 ± 0.2</td></tr><tr><td>Scale height aspect ratio</td><td>h<sub>0</sub></td><td>...</td><td>0.0001–0.5</td><td>0.079 ± 0.005</td><td>0.084 ± 0.006</td></tr><tr><td>HG asymmetry parameter</td><td>g</td><td>...</td><td>0–1</td><td>0.285 ± 0.012</td><td>0.65 ± 0.03</td></tr><tr><td>Line-of-sight inclination</td><td>φ</td><td>°</td><td>45–85</td><td>72.2 ± 0.4</td><td>68.3<sup>+0.7</sup><sub>-0.8</sub></td></tr><tr><td>Position angle</td><td>θ<sub>PA</sub></td><td>°</td><td>105°–145°</td><td>127.0 ± 0.3</td><td>127.1 ± 0.3</td></tr><tr><td>Sky plane offset</td><td>ΔX<sub>1</sub></td><td>au</td><td>−R<sub>1</sub>–R<sub>1</sub></td><td>1.6 ± 0.6</td><td>1.4 ± 0.6</td></tr><tr><td>Line-of-sight offset</td><td>ΔX<sub>2</sub></td><td>au</td><td>−R<sub>1</sub>–R<sub>1</sub></td><td>2.1 ± 1.6</td><td>−5.7<sup>+2.1</sup><sub>-2.2</sub></td></tr><tr><td>Chi-Squared</td><td>χ<sup>2</sup></td><td>...</td><td>...</td><td>2648</td><td>2670</td></tr><tr><td>Reduced Chi-Squared</td><td>χ<sup>2</sup><sub>red</sub></td><td>...</td><td>...</td><td>0.89</td><td>0.90</td></tr></tbody></table>

is not included in the table. Plots of the marginalized and joint posterior distributions derived from the MCMC chains can be seen in appendix Figure 7 and appendix Figure 8, for the HG only and HG + Rayleigh models, respectively. The posteriors for each parameter are all single-peaked and approximately Gaussian, with the exception of the outer radius, $R_2$ that displays a slightly elongated tail toward higher values. In the HG model, slight degeneracies are found between $\phi$ and $\Delta X_1$, and $\phi$ and $\Delta X_2$. The HG + Rayleigh model displays addition small degeneracies between $g$ and $\phi$, as well as $g$ and $\Delta X_2$. However, all parameters appear to be well constrained in the marginalized distributions. Images of the best-fit models for both the HG fit and the HG + Rayleigh fit can be seen in Figure 6, alongside an image of the residuals of each model subtracted from the $Q_r$ image.

Table 3 displays the root mean square (rms) pixel values, as well as the fifth and ninety-fifth percentile pixel values (as a proxy for dynamic range) for the $Q_r$ image, the $U_r$ image, and the two residual images. The measurements have been split into an inner region (inside 0"/275) and outer region (between 0"/3 and 1"/3). The inner region is dominated by residual systematics and has very little disk flux. As a result, subtracting the models from the $Q_r$ image has little effect on this region, as demonstrated by the fact that the percentiles and the rms in the residual images remain close to the original $Q_r$ values. The differences between the noise in $Q_r$ images and the $U_r$ images at these separations may be due to how the residual instrumental polarization affects each image. However, this remains to be confirmed. In the outer region, after the models have been subtracted all three quantities converge to the $U_r$ values, indicating that we are successfully subtracting off the disk flux with the residuals consistent with pure noise (i.e., the $U_r$ image).

The best-fit HG (HG + Rayleigh) model reveals a disk inclination of 72°.2 ± 0.4 (68°.3+0.7-0.8), with a relatively steep radial density power-law index of 2.25 ± 0.15 (2.2 ± 0.2), that extends from an inner radius of 81 ± 2 au (78.9 ± 0.8 au) to an outer radius of 216 ± 16 au (211+21-15 au). The best-fit inclinations are roughly consistent with that determined by measuring the aspect ratio in the STIS disk image (~75°; D. Padgett et al. 2016, in preparation). The outer radius is constrained by the level of forward scattering seen in the $Q_r$ image toward the NE, which quickly reaches the noise floor. Therefore, this fit parameter is governed by the sensitivity of these observations, and represents a lower limit on the true disk outer radius.

<table><caption>**Table 3**
Pixel Statistics for the Data and Model Residual Images</caption><thead><tr><th rowspan="2">Frame</th><th colspan="3">Inside 0"/275</th><th colspan="3">Between 0"/275 and 1"/3</th></tr><tr><th>5%-tile</th><th>95%-tile</th><th>rms</th><th>5%-tile</th><th>95%-tile</th><th>rms</th></tr></thead><tbody><tr><td>Q<sub>r</sub></td><td>-26.0</td><td>18.3</td><td>13.6</td><td>-3.1</td><td>8.6</td><td>3.8</td></tr><tr><td>U<sub>r</sub></td><td>-18.5</td><td>24.1</td><td>12.8</td><td>-3.9</td><td>3.4</td><td>2.4</td></tr><tr><td>HG Residuals</td><td>-27.4</td><td>14.4</td><td>13.4</td><td>-4.1</td><td>3.3</td><td>2.4</td></tr><tr><td>HG+Rayleigh Residuals</td><td>-26.2</td><td>15.6</td><td>13.2</td><td>-3.8</td><td>3.7</td><td>2.4</td></tr></tbody></table>

Note. All values in this table are presented in raw data units (i.e., ADU/coadd).

Indeed, the disk is seen to extend to a radius greater than 500 au in the HST/STIS images. However, it is possible that the HST/STIS observations probe a different dust population.

The results of the disk fitting also indicate that the disk is offset from the star along the projected semimajor axis, $\Delta X_1$, by 1.6 ± 0.6 au (1.4 ± 0.6 au). This feature is consistent between both models, with both fits finding an offset with a significance slightly below 3σ. The direction of this offset is consistent with the direction of the brightness asymmetries seen in Figure 1 and causes a faint brightness asymmetry in both model images in Figure 6.

Visually, there are two major differences between the two models. First, there is an apparent deficit of light near the projected semimajor axis in the HG + Rayleigh model, which is due to the low polarization fraction at small (and large) scattering angles of the Rayleigh function. The deficit is in an area that coincides with a region of reduced S/N in Figure 2 and it is therefore likely that this area was down-weighted in the MCMC fitting relative to the ansae. In Figure 1, the brightness appears to increase toward the NE, with a maximum near the semiminor axis and no apparent deficit. This may be a real feature of the disk, but could also be possibly due to residual systematics in the image. Indeed, there are strong residuals at similar separations along the SW semiminor axis for both models. Comparing the residuals between the two models in Figure 6, there appears to be additional flux at the location of the deficit in the HG + Rayleigh model. However, the nearly identical reduced-$\chi^2$ value indicates that both models fit the data equally well.

The second major difference between the two models is the brightness of the backside of the disk, which is barely visible in
the HG + Rayleigh model but appears in the image of the HG model. A close examination of the HG + Rayleigh model residual image, reveals a faint brightness at the location of the backside in the HG model that are not apparent in the HG-only model residuals. These residuals may indicate a very low S/N detection of the backside of the disk in the *Q*<sub>r</sub> image.

Between the two models, the $\Delta X_2$ offset varies in both magnitude and sign. In the HG-only best-fit model, $\Delta X_2$ is positive, indicating that the star is offset away from the observer along the inclination vector, relative to the disk center. The HG + Rayleigh model is best fit by a negative offset, where the star is offset toward the observer. This discrepancy illuminates a degeneracy between this offset and the exact form of the polarized scattering function. We therefore consider this parameter to remain poorly constrained. Future studies that are able to constrain the grain scattering properties or image the backside of the disk will be able to elucidate this remaining unknown.

## 5. DISCUSSION

The dust seen in scattered-light images of debris disks is thought to have originated in disks or belts of planetesimals, where collisional cascades grind kilometer-sized bodies down to micron-sized dust (see, Wyatt 2008, and reference therein). In order to initiate these cascades, the constituent planetesimals must be dynamically stirred such that their eccentricities reach a high enough level (on the order of $10^{-3}-10^{-2}$) to allow their orbits to cross and for their collisions to be destructive. Stirring mechanisms include: self-stirring, where objects on the order of 2000 km located inside the belts induce the cascade from within (e.g., Kenyon & Bromley 2004); planet-induced stirring, where a nearby planet excites the disk (e.g., Mustill & Wyatt 2009); or dynamical interaction with a passing star (e.g., Kenyon & Bromley 2002). With the exception of the stellar fly-by scenario, the strength of all of these mechanisms should diminish with stellar age. As a given system reaches a steady-state configuration it cools dynamically and the collision rate slows. Therefore, the scattered-light luminosity should dim with age as the small grains are removed via radiation pressure, Poynting-Robertson drag and/or stellar winds, and can no longer be replenished through collisions. Note that for disks of all ages, large collisions or other transient events may cause a temporary increase in dust production and create a short-term increase in disk brightness (e.g., Wyatt et al. 2007).

The only published age estimates of HD 157587 are presented by Feltzing et al. (2001) and Casagrande et al. (2011). Feltzing et al. (2001) fit metallicity, effective temperature, and absolute magnitude to the outputs of a rapid stellar evolution algorithm (Hurley et al. 2000) based on the evolutionary tracks produced by Pols et al. (1998). They estimate HD 157587 to have an age of $2.2 \pm 0.5$ Gyr. Casagrande et al. (2011) use a Bayesian analysis to fit effective temperature, metallicity, and Johnson V magnitude to the Padova evolutionary tracks (Bertelli et al. 2008, 2009) and the BASTI isochrones (Pietrinferni et al. 2004, 2006, 2009). They find an age of $3.0_{-1.5}^{+1.7}$ Gyr and $3.0_{-1.4}^{+1.0}$ Gyr at 95% confidence for the Padova and BASTI models, respectively.

Such advanced ages suggest that the star has evolved along the HR diagram away from the zero age main sequence (ZAMS). Alternatively, this offset from the ZAMS may be indicative of a pre-main-sequence star moving toward the ZAMS, rather than away from it. However, considering the

timescales of evolution along the pre-main-sequence tracks (~12 Myr for an F5 star, Siess et al. 2000) compared to main-sequence tracks (on the order of gigayears), it is more likely that the star has been found as it moves away from the ZAMS.

Conversely, applying the star's proper motion, radial velocity, and parallax to the BANYAN II webtool (Gagné et al. 2014; Malo et al. 2013)<fnref n="32" /> indicates a 91% probability that the star is a young (<1 Gyr) field star, and a 9% probability that the star is an old field star. Additionally, with velocities of $[U, V, W] = [-7, -17, -8]$ km s⁻¹ (Holmberg et al. 2009), we find that HD 157587 could be kinematically associated with several relatively young moving groups. For example, in Chereul et al. (1999) the Pleiades stream 2–5 has $[U, V, W] = [-12.0 \pm 5.3, -21.6 \pm 4.7, -5.3 \pm 5.9]$ km s⁻¹, whereas the Centaurus–Lupus stream 2–12 has $[U, V, W] = [-12.4 \pm 6.1, -16.5 \pm 4.6, -7.4 \pm 3.1]$ km s⁻¹. In Asiain et al. (1999), the Pleiades moving group B2 has $[U, V, W] = [-10.7 \pm 5.3, -18.8 \pm 3.7, -5.6 \pm 2.2]$ km s⁻¹, whereas members of Lower Centaurus Crux have $[U, V, W] = [-6.8 \pm 4.7, -18.5 \pm 6.5, -6.4 \pm 1.7]$ km s⁻¹. In all of these cases, the age distributions determined by various methods lie in the range of $10^7-10^9$ years. Thus, the stellar kinematics possibly support a relatively young age for HD 157587, while its photometry supports a much older age.

If HD 157587 is truly 2.5–3 Gyr old, then its disk's optical/NIR emission would be unusually bright when compared to the overall population of imaged debris disks (for a good summary, see Figure 1 from Choquet et al. 2016). In fact, the disk around HD 157587 would represent one of only five debris disks imaged in scattered light with ages greater than one gigayear. The other four old disks seen in scattered light are HD 207129 (Krist et al. 2010), HD 202629 (Krist et al. 2012), HD 53143 (Kalas et al. 2006), and HD 10647 (q¹ Eri; Stapelfeldt et al. 2007). With a spectral type of F5V, HD 157587 is the earliest type star of the five, making it potentially the earliest type star older than 1 Gyr with a debris disk seen in scattered light.

Of the other four old scattered-light debris disks (i.e., >1 Gyr), HD 10647 is the most similar to HD 157587. HD 10647 is an F8V star with a debris disk at a radial distance of ~85 au (Liseau et al. 2010). The disk has been imaged in scattered light with HST/ACS (Stapelfeldt et al. 2007) and in the infrared with Herschel/PACS (Liseau et al. 2010). As seen here for HD 157587, the HST images show evidence of a brightness asymmetry between the two ansae of the disk. Interestingly, HD 10647 is known to host a Jupiter mass planet at 2 au. However, this planet is at too great a distance from the disk to have significant dynamical influence on the disk. Liseau et al. (2008) posit that the disk asymmetry suggests the presence of a second planet at larger distances.

In our images of the HD 157587 disk, we see a slight brightness asymmetry between the SE and NW ansae of the disk that we have modeled as being caused by a stellocentric offset. If the offset is due to perturbations from a substellar companion, it is highly unlikely that it is one of the three point sources imaged in total intensity. Using the COND evolutionary models (Baraffe et al. 2003) and assuming an age of 3 Gyr, we find their luminosities (flux ratios) correspond to brown dwarf masses between 30 and 40 $M_{Jup}$. If these three objects reside within the plane of the disk, their deprojected

<footnote n="32">³² http://www.astro.umontreal.ca/~gagne/banyanII.php</footnote>
separations correspond to stellocentric radii between 210 and 380 au, assuming circular orbits. Considering the smooth radial extent of the disk seen in the STIS images, we deem such an alignment unlikely, as the presence of such massive companions would cause significant disturbances to the disk morphology. Though it is possible that the expulsion timescale for the STIS grains are short enough to hide any such feature.

Future imaging observations that reveal the relative motion of these three sources are required to understand their true relationship to HD 157587, if any. However, it is worth noting that the likelihood of a field star flyby is rare—less than 1% in 100 Myr of a 500 au approach (Kenyon & Bromley 2002)—and considering HD 157587's proximity to the galactic plane, these three objects are likely background sources. If HD 157587 is a younger star, as suggested by the stellar kinematics, it is possible that a lower mass perturber at smaller angular separations is still bright enough in thermal emission to be detected via direct imaging.

## 6. CONCLUSION

Using GPI, we have imaged the dust ring around HD 157587 in H-band polarized intensity. The image reveals an inclined disk that appears to be cleared of material inside of a projected major axis of ~80 au. The FOV of our observations overlaps with the inner regions of previous STIS images of the disk, and our analysis returns a similar disk inclination to that derived with the STIS data. The disk has a strong polarized brightness asymmetry in the NE-SW direction, where we interpreted the bright side of the disk to be tilted toward the observer. A similar brightness asymmetry has been seen in polarized observations of a number of other recently imaged disks, suggestive of similar grain compositions, size distributions, and/or dust grain morphologies. Future detailed studies of these disks' dust composition that include multicolor observations or polarization fraction measurements will be able to further explore the similarities and differences of their dust grain populations.

A second, weaker, brightness asymmetry is seen between the two ansae that could be due to a stellocentric offset in the plane of the sky. To test this hypothesis we used Bayesian MCMC methods to fit the polarized disk image to two disk models, one that used an HG polarized scattering phase function and one that combined an HG function with a Rayleigh scattering phase function. Both models reveal an offset dust disk with an inner radius of 80 au and an inclination of about 70°. The center of the disk is found to be offset approximately 1.5 au from the star's location in the plane of the sky and both models reproduce the brightness asymmetry between the two ansae. This offset could be confirmed with longer wavelength imaging using ALMA, which would trace thermal emission and therefore have less of a dependence on the scattering properties of the grains.

In general, the two model fits return similar disk properties, with the exception of $\Delta X_2$, the offset in the disk plane. We find that the form of the polarized scattered phase function is degenerate with the magnitude and direction of this offset and

without further information on the form of the scattering phase function this value will remain poorly constrained.

The total intensity observations are dominated by stellar residuals at the location of the disk in the polarized intensity image and no disk was recovered. However, three point sources were recovered. Considering HD 157587's proximity to the galactic plane and the positions of the point sources relative to the disk, we consider these point sources to be background objects. Nonetheless, follow-up observations are required to confirm this proposition.

The currently published ages of the system that rely on stellar evolutionary tracts indicate an age well over 1 Gyr old. However, such an evolved age is at odds with the stellar kinematics. The stellocentric offset suggests that this system has a complicated dynamical history and may harbour one or more unseen planets. This notion is reinforced by the similarities between HD 157587's stellar properties and disk morphology, and those of the RV planet host HD 10647. If the stellocentric offset is due to perturbations by one or more planets, further detailed study of the system's debris disk will be required to thoroughly characterize the system; the system's advanced age would make it ill-suited for direct imaging planet searches and radial velocity measurements would require prohibitively long time baselines. On the other hand, if the disk is younger than 1 Gyr, as implied by the stellar kinematics, then it presents itself as a prime target for deeper direct imaging observations, which may be able to image the disk's perturber.

This work was based on observations obtained at the Gemini Observatory, which is operated by the Association of Universities for Research in Astronomy, Inc., under a cooperative agreement with the National Science Foundation (NSF) on behalf of the Gemini partnership: the NSF (United States), the National Research Council (Canada), CONICYT (Chile), the Australian Research Council (Australia), Ministério da Ciência, Tecnologia e Inovação (Brazil) and Ministerio de Ciencia, Tecnología e Innovación Productiva (Argentina).

This work was supported in part by NASA's NEXSS program, grant number NNX15AD95G. Portions of this work were performed under the auspices of the U.S. Department of Energy by Lawrence Livermore National Laboratory under Contract DE-AC52-07NA27344. Portions of this work were also carried out with the support of NSF grants AST-1413718 and AST-1518332, and NASA grant NNX15AC89G. This research used resources of the National Energy Research Scientific Computing Center, a DOE Office of Science User Facility supported by the Office of Science of the U.S. Department of Energy under Contract No. DE-AC02-05CH11231.

Facility: Gemini:South(GPI).

## APPENDIX

The marginalized and joint probability distributions of the disk model parameters from the MCMC fitting of the two disk models from Section 4 can be seen in Figures 7 and 8.
<figure><img src="image_9.png" /><figcaption>**Figure 7.** Posterior distributions of the model parameters from MCMC disk model fitting to the *Q* disk image with the HG only model. The diagonal histograms show the posterior distributions of each parameter marginalized across all the other parameters. In each plot, the dashed lines indicate the 16%, 50%, and 84% percentiles. The off-diagonal plots display the joint probability distributions with contour levels at the same percentiles. The normalization term, *N*₀, has been excluded from this plot and is considered a nuisance parameter.</figcaption></figure>
<figure><img src="image_10.png" /><figcaption>**Figure 8.** Posterior distributions of the model parameters from MCMC disk model fitting to the *Q*<sub>r</sub> disk image with the HG + Rayleigh model. The diagonal histograms show the posterior distributions of each parameter marginalized across all the other parameters. In each plot, the dashed lines indicate the 16%, 50%, and 84% percentiles. The off-diagonal plots display the joint probability distributions with contour levels at the same percentiles. The normalization term, *N*<sub>0</sub>, has been excluded from this plot and is considered a nuisance parameter.</figcaption></figure>

REFERENCES

Asiain, R., Figueras, F., & Torra J. 1999, A&A, 350, 434

Baraffe, I., Chabrier, G., Barman, T. S., Allard, F., & Hauschildt, P. H. 2003, A&A, 402, 701

Bertelli, G., Girardi, L., Marigo, P., & Nasi, E. 2008, A&A, 484, 815

Bertelli, G., Nasi, E., Girardi, L., & Marigo, P. 2009, A&A, 508, 355

Burrows, C. J., Krist, J. E., Stapelfeldt, K. R. & WFPC2 Investigation Definition Team 1995, BAAS, 27, 1329

Casagrande, L., Schönrich, R., Asplund, M., et al. 2011, A&A, 530, A138

Chereul, E., Crézé, M., & Bienaymé, O. 1999, A&AS, 135, 5

Choquet, É., Perrin, M. D., Chen, C. H., et al. 2016, ApJL, 817, L2

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

Esposito, T. M., Fitzgerald, M. P., Graham, J. R., et al. 2016, AJ, 152, 85

Feltzing, S., Holmberg, J., & Hurley, J. R. 2001, A&A, 377, 911

Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306

Gagné, J., Lafrenière, D., Doyon, R., Malo, L., & Artigau, É. 2014, ApJ, 783, 121

Graham, J. R., Kalas, P. G., & Matthews, B. C. 2007, ApJ, 654, 595

Hedman, M. M., & Stark, C. C. 2015, ApJ, 811, 67

Holmberg, J., Nordström, B., & Andersen, J. 2009, A&A, 501, 941

Hung, L.-W., Duchêne, G., Arriaga, P., et al. 2015, ApJL, 815, L14

Hurley, J. R., Pols, O. R., & Tout, C. A. 2000, MNRAS, 315, 543

Kalas, P., Duchene, G., Fitzgerald, M. P., & Graham, J. R. 2007, ApJL, 671, L161

Kalas, P., Graham, J. R., & Clampin, M. 2005, Natur, 435, 1067

Kalas, P., Graham, J. R., Clampin, M. C., & Fitzgerald, M. P. 2006, ApJL, 637, L57

Kalas, P. G., Rajan, A., Wang, J. J., et al. 2015, ApJ, 814, 32
Kenyon, S. J., & Bromley, B. C. 2002, AJ, 123, 1757

Kenyon, S. J., & Bromley, B. C. 2004, AJ, 127, 513

Krist, J. E., Stapelfeldt, K. R., Bryden, G., et al. 2010, AJ, 140, 1051

Krist, J. E., Stapelfeldt, K. R., Bryden, G., & Plavchan, P. 2012, AJ, 144, 45

Liseau, R., Eiroa, C., Fedele, D., et al. 2010, A&A, 518, L132

Liseau, R., Risacher, C., Brandeker, A., et al. 2008, A&A, 480, L47

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

Maire, J., Perrin, M. D., Doyon, R., et al. 2012, Proc. SPIE, 8451, 3

Malo, L., Doyon, R., Lafrenière, D., et al. 2013, ApJ, 762, 88

Mawet, D., Milli, J., Wahhaj, Z., et al. 2014, ApJ, 792, 97

McDonald, I., Zijlstra, A. A., & Boyer, M. L. 2012, MNRAS, 427, 343

Millar-Blanchaer, M. A., Graham, J. R., Pueyo, L., et al. 2015, ApJ, 811, 18

Mouillet, D., Larwood, J. D., Papaloizou, J. C. B., & Lagrange, A. M. 1997, MNRAS, 292, 896

Mustill, A. J., & Wyatt, M. C. 2009, MNRAS, 399, 1403

Padgett, D., & Stapelfeldt, K. 2015, in IAU Symp. 314, Young Stars and Planets Near the Sun, ed. J. H. Kastner, B. Stelzer, & S. A. Metchev, 175

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, 3

Pietrinferni, A., Cassisi, S., Salaris, M., & Castelli, F. 2004, ApJ, 612, 168

Pietrinferni, A., Cassisi, S., Salaris, M., & Castelli, F. 2006, ApJ, 642, 797

Pietrinferni, A., Cassisi, S., Salaris, M., Percival, S., & Ferguson, J. W. 2009, ApJ, 697, 275

Pols, O. R., Schröder, K.-P., Hurley, J. R., Tout, C. A., & Eggleton, P. P. 1998, MNRAS, 298, 525

Quillen, A. C. 2006, MNRAS, 372, L14

Schmid, H. M., Joos, F., & Tschan, D. 2006, A&A, 452, 657

Siess, L., Dufour, E., & Forestini, M. 2000, A&A, 358, 593

Soummer, R., Pueyo, L., & Larkin, J. 2012, ApJL, 755, L28

Stapelfeldt, K., Krist, J., Bryden, G., & Chen, C. 2007, in *In the Spirit of Bernard Lyot: The Direct Detection of Planets and Circumstellar Disks in the 21st Century*, ed. P. Kalas, 47

Telesco, C. M., Fisher, R. S., Piña, R. K., et al. 2000, ApJ, 530, 329

van Leeuwen, F. 2007, A&A, 474, 653

Wang, J. J., Rajan, A., Graham, J. R., et al. 2014, Proc. SPIE, 9147, 55

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

Wyatt, M. C. 2008, ARA&A, 46, 339

Wyatt, M. C., Dermott, S. F., Telesco, C. M., et al. 1999, ApJ, 527, 918

Wyatt, M. C., Smith, R., Greaves, J. S., et al. 2007, ApJ, 658, 569