Exploring the Brown Dwarf Desert: New Substellar Companions from the SDSS-III MARVELS Survey

TL;DR

This study reports 10 new brown dwarf companions and 2 low-mass stellar companions around solar-type stars, providing new data to understand the brown dwarf desert phenomenon and its implications for substellar object formation.

Contribution

The paper introduces new substellar companions discovered via the MARVELS survey, utilizing an improved data pipeline to reduce errors and enhance detection of brown dwarfs.

Findings

10 brown dwarf companions with masses 13-76 M_Jup within 1 AU

Estimated brown dwarf occurrence rate of 0.56% around solar-type stars

Improved detection pipeline reduces systematic errors in radial velocity data

Abstract

Planet searches using the radial velocity technique show a paucity of companions to solar-type stars within ~5 AU in the mass range of ~10 - 80 M$_{\text{Jup}}$. This deficit, known as the brown dwarf desert, currently has no conclusive explanation. New substellar companions in this region help asses the reality of the desert and provide insight to the formation and evolution of these objects. Here we present 10 new brown dwarf and two low-mass stellar companion candidates around solar-type stars from the Multi-object APO Radial-Velocity Exoplanet Large-Area Survey (MARVELS) of the Sloan Digital Sky Survey III (SDSS-III). These companions were selected from processed MARVELS data using the latest University of Florida Two Dimensional (UF2D) pipeline, which shows significant improvement and reduction of systematic errors over previous pipelines. The 10 brown dwarf companions range in mass from ~13 to 76 M$_{\text{Jup}}$ and have orbital radii of less than 1 AU. The two stellar companions have minimum masses of ~98 and 100 M$_{\text{Jup}}$. The host stars of the MARVELS brown dwarf sample have a mean metallicity of [Fe/H] = 0.03 $\pm$ 0.08 dex. Given our stellar sample we estimate the brown dwarf occurrence rate around solar-type stars with periods less than ~300 days to be ~0.56%.

Tables

Table 1. Table 1: MARVELS Host Star Parameters

Starname	T${}_{\text{eff}}$	[Fe/H]	log g	M⋆	R⋆	V${}_{\text{mag}}$
	(K)	(dex)	(cgs)	(M⊙)	(R⊙)	(GSC 2.3)
TYC 0123-00549-1† (MCKGS1-50)	5540 $\pm$ 39	0.20 $\pm$ 0.07	4.48 $\pm$ 0.23	1.02 $\pm$ 0.08	0.95 $\pm$ 0.27	10.3
TYC 0173-02410-1	5470 $\pm$ 112	0.18 $\pm$ 0.09	4.25 $\pm$ 0.18	1.06 $\pm$ 0.10	1.26 $\pm$ 0.30	10.5
TYC 4955-00369-1	6190 $\pm$ 143	-0.43 $\pm$ 0.09	4.55 $\pm$ 0.24	1.03 $\pm$ 0.09	0.90 $\pm$ 0.26	10.2
TYC 1451-01054-1	5980 $\pm$ 140	-0.34 $\pm$ 0.09	4.60 $\pm$ 0.24	0.99 $\pm$ 0.08	0.83 $\pm$ 0.24	9.9
TYC 3469-00492-1	5285 $\pm$ 82	0.36 $\pm$ 0.06	3.78 $\pm$ 0.13	1.29 $\pm$ 0.12	2.39 $\pm$ 0.46	11.3
GSC 03467-00030† (MCKGS1-135)	5525 $\pm$ 71	-0.17 $\pm$ 0.08	4.39 $\pm$ 0.23	0.96 $\pm$ 0.09	1.03 $\pm$ 0.30	12.2
TYC 3091-00716-1† (MCKGS1-153)	5614 $\pm$ 29	-0.09 $\pm$ 0.06	4.61 $\pm$ 0.09	0.95 $\pm$ 0.06	0.80 $\pm$ 0.09	10.0
TYC 3097-00398-1	5340 $\pm$ 98	0.16 $\pm$ 0.07	3.85 $\pm$ 0.15	1.21 $\pm$ 0.12	2.14 $\pm$ 0.47	11.3
TYC 3130-00160-1† (MCUF1-11)	5315 $\pm$ 44	0.33 $\pm$ 0.06	4.23 $\pm$ 0.19	1.06 $\pm$ 0.10	1.28 $\pm$ 0.33	11.0
TYC 3547-01007-1	6090 $\pm$ 190	0.06 $\pm$ 0.11	4.57 $\pm$ 0.27	1.12 $\pm$ 0.10	0.91 $\pm$ 0.30	10.5
TYC 3556-03602-1	6295 $\pm$ 200	0.11 $\pm$ 0.09	4.35 $\pm$ 0.30	1.25 $\pm$ 0.14	1.23 $\pm$ 0.48	10.6
TYC 3148-02071-1	6270 $\pm$ 271	0.09 $\pm$ 0.14	4.13 $\pm$ 0.27	1.33 $\pm$ 0.18	1.63 $\pm$ 0.62	9.5

Table 2. Table 2: MARVELS Companion Parameters

Star	MARVELS	M sini	Period	Ecc	K
	designation	(M${}_{\text{Jup}}$)	(days)		(m s-1)
TYC 0123-00549-1 (MCKGS1-50)	MARVELS-8b	37.7 $\pm$ 0.4	36.045 $\pm$ 0.005	0.101 $\pm$ 0.004	2241 $\pm$ 8
TYC 0173-02410-1	MARVELS-9b	76.0 $\pm$ 5.1	216.5 $\pm$ 1.8	0.37 $\pm$ 0.05	2537 $\pm$ 210
TYC 4955-00369-1	MARVELS-10b	14.8 $\pm$ 2.2	217.3 $\pm$ 1.1	0.53 $\pm$ 0.15	574 $\pm$ 128
TYC 1451-01054-1	MARVELS-11b	52.1 $\pm$ 0.6	11.6121 $\pm$ 0.0003	0.000 $\pm$ 0.002	4551 $\pm$ 12
TYC 3469-00492-1	MARVELS-12b	32.8 $\pm$ 1.1	252.5 $\pm$ 0.1	0.60 $\pm$ 0.01	1093 $\pm$ 35
GSC 03467-00030 (MCKGS1-135)	MARVELS-13b	41.8 $\pm$ 2.9	147.6 $\pm$ 0.3	0.50 $\pm$ 0.04	1868 $\pm$ 167
TYC 3091-00716-1 (MCKGS1-153)	MARVELS-14b	100.2 $\pm$ 1.4	269.1 $\pm$ 0.3	0.417 $\pm$ 0.003	3359 $\pm$ 14
TYC 3097-00398-1	MARVELS-15b	97.6 $\pm$ 1.8	95.54 $\pm$ 0.04	0.223 $\pm$ 0.004	3730 $\pm$ 17
TYC 3130-00160-1 (MCUF1-11)	MARVELS-16b	42.7 $\pm$ 0.7	31.656 $\pm$ 0.007	0.690 $\pm$ 0.006	3552 $\pm$ 42
TYC 3547-01007-1	MARVELS-17b	13.1 $\pm$ 1.6	103.4 $\pm$ 0.9	0.24 $\pm$ 0.08	537 $\pm$ 64
TYC 3556-03602-1	MARVELS-18b	24.4 $\pm$ 2.4	117.6 $\pm$ 1.1	0.46 $\pm$ 0.09	968 $\pm$ 83
TYC 3148-02071-1	MARVELS-19b	17.3 $\pm$ 2.4	122.2 $\pm$ 1.7	0.53 $\pm$ 0.10	683 $\pm$ 82

Table 3. Table 3: Vetting Criteria Dependencies

Criteria	original	minus 1$\sigma$			plus 1$\sigma$
	dependence	dependence	# pass	#BD pass	dependence	# pass	#BD pass
1	5%	3.30%	84	27	6.71%	84	27
2	$e>0.934$	$e>0.62$	76	23	$e>0.99$	86	27
3			84	27		79	25
4	$3+3\left(1-V_N\right)e$	$1.98+1.98\left(1-V_N\right)e$	85	27	$4.02+4.02\left(1-V_N\right)e$	80	27
5	$3+3\left(1-V_N\right)e$	$1.98+1.98\left(1-V_N\right)e$	84	27	$4.02+4.02\left(1-V_N\right)e$	84	27
6	$3+3\left(1-V_N\right)e$	$1.98+1.98\left(1-V_N\right)e$	95	31	$4.02+4.02\left(1-V_N\right)e$	77	24
7	$3+3\left(1-V_N\right)e$	$1.98+1.98\left(1-V_N\right)e$	84	27	$4.02+4.02\left(1-V_N\right)e$	80	24
8			43	14		85	28
9	$e$ $>$ 0.4	$e$ $>$ 0.26	78	26	$e$ $>$ 0.54	86	28
10	$\Delta T\ge 2P$	$\Delta T\ge 1.32P$	94	33	$\Delta T\ge 2.68P$	76	21
11	$\ge$ 1/3	$\ge$ 0.22	89	32	$\ge$ 0.45	76	24

Equations

$$\mu =\sqrt{(cosd\ast {\mu }_r{)}^2+{\mu }_d^2}$$

$${\text{RPM}}_J=J+5log(\mu ),$$

$$\begin{array}{r}y=-58+313.42(J-H)-583.6(J-H{)}^2\\ +473.18(J-H^3-141.25(J-H{)}^4,\end{array}$$

$$\Sigma_{RV}=\text{stddev}\Bigg{(}\frac{\textbf{v}-\widetilde{\text{v}}}{\sigma_{\textbf{v}}}\Bigg{)}\geq 2.5,$$

$$U_{N}=\frac{N}{N-1}\bigg{[}1-\sum_{i=1}^{N}(\phi_{i+1}-\phi_{i})^{2}\bigg{]},$$

$$V_{N}=\frac{N}{N-1}\bigg{[}1-\sum_{i=1}^{N}(\nu_{i+1}-\nu_{i})^{2}\bigg{]},$$

$${\chi }_{\text{mod}}^2=\frac{{\chi }_{\text{red}}^2}{\sqrt{U_N{V}_N}},$$

Full text

††pubyear: 2017††pagerange: Exploring the Brown Dwarf Desert: New Substellar Companions from the SDSS-III MARVELS Survey–B

1 Introduction

Brown dwarfs (BDs) are substellar objects ranging in mass between low-mass main-sequence stars and giant planets. They are generally defined as having masses above $\sim$13 M${}_{\text{Jup}}$, the minimum mass required for an object to ignite deuterium burning in its core, and below $\sim$80 M${}_{\text{Jup}}$, or 0.075 M⊙, above which a dense object becomes sufficiently massive to fuse atomic hydrogen into helium in its core (Boss et al., 2003). The existence of BDs was first suggested by Kumar (1962), on the basis of theoretical calculations, and the realization that as stars of sufficiently low mass condense out of the interstellar medium, they will become dense enough to be supported by electron degeneracy pressure before they become hot enough to fuse hydrogen in their core. However, there was not a confirmed observation of a BD until 1995 (Rebolo et al., 1995; Nakajima et al., 1995). Technological advancements continue to produce new BD discoveries with increasing detail. These new data give a statistical sense of these objects and place observational constraints on their theoretical properties and proposed formation mechanisms.

One of the most important observational constraints on BD formation is known as the “brown dwarf desert”, or the lack of companions within the mass range of $\sim$10-80 M${}_{\text{Jup}}$ with orbital radii less than $\sim$5 AU around solar-type stars (Marcy & Butler, 2000). Grether & Lineweaver (2006) investigated the significance of the brown dwarf desert by analyzing companions with orbital periods less than 5 years around Sun-like stars. Grether & Lineweaver (2006) found that $\sim$16% of Sun-like stars have close (P $<$ 5 yr) companions more massive than Jupiter, but $<$1% are BDs. The minimum number of companions per unit interval in log mass, or the driest part of the desert, is at M $\sim$ 31 M${}_{\text{Jup}}$ (Grether & Lineweaver, 2006). Ma & Ge (2014) created a catalogue of all known brown dwarf companions around solar-type stars, which totaled 62 BDs. Through statistical analysis Ma & Ge (2014) identified a depletion of brown dwarfs with masses between 35 and 55 M${}_{\text{Jup}}$ and periods $<$ 100 days, representing the ‘driest’ part in the brown dwarf desert.

Troup et al. (2016) consider stellar and substellar companion candidates in the Sloan Digital Sky Survey III (SDSS-III; Eisenstein et al. (2011)) Apache Point Observatory Galactic Evolution Experiment (APOGEE-1; Majewski et al. (2015)), which observed 146,000 stars with typical radial velocity (RV) precision of $\sim$100-200 m s*-1*. Of these, 14,840 had at least eight visits; this subset had a median baseline of 384 days, median of 13 visits, and a median Signal-to-Noise ratio (S/N) of 12.2 per visit. With a final count of 382 candidate companions (including 112 BD companion candidates, 71 orbiting evolved stars), Troup et al. (2016) find evidence for the brown dwarf desert, but only out to $\sim$0.5-1 AU, with no indication of a desert beyond. However, of the 382 stars in their sample, only 36 are solar-type (main-sequence with 5000K $<$ T${}_{\text{eff}}$ $<$ 6000K) stars. Troup et al. (2016) argue that this result agrees with the work of Guillot et al. (2014) who find that the the BD desert may not be as “dry" as initially thought when considering stars more massive than the Sun, as well as evidence that the BD desert is a special case for solar-mass stars (Duchêne & Kraus, 2013). Duchêne & Kraus (2013) suggest tidal interaction with host stars may shape the brown dwarf desert. While assessing the characteristic timescale for orbital decay, Damiani & Díaz (2016) find that F-type stars may host massive companions for a significantly longer time than G-type stars for orbital periods less than five days; however, they conclude that brown dwarf occurrence should largely be unaffected by tidal decay for orbital periods greater than ten days, regardless of host mass.

Debate continues on valid formation mechanisms of BDs and whether or not different formation mechanisms dominate in different mass, period, and eccentricity regimes. Three main mechanisms are advocated in the field: (i) the ejection of protostellar embryos from their natal cloud cores (Boss, 2001; Reipurth & Clarke, 2001), (ii) the fragmentation of protostellar discs (Stamatellos & Whitworth, 2009; Stamatellos et al., 2011), or (iii) collapsing molecular cloud cores (Hennebelle & Chabrier, 2013; Hopkins, 2012), i.e., the same process as Sun-like stars.

Through their statistical study Ma & Ge (2014) reported that BD companions with minimum masses below 42.5 M${}_{\text{Jup}}$ have an eccentricity distribution consistent with that of massive planets, while those with masses above 42.5 M${}_{\text{Jup}}$ have an eccentricity distribution consistent with that of binaries. Ma & Ge (2014) also find that host stars of BD companions are not metal rich and have a significantly different metallicity distribution compared to host stars of giant planets, suggesting a formation mechanism different from the core-accretion (Ida & Lin, 2004) scenario. Ma & Ge (2014) conclude that BD companions with minimum masses below 42.5 M${}_{\text{Jup}}$ likely form in a protoplanetary disc through the disc gravitational instability scenario, and their eccentricity is excited through scattering, while BD companions with minimum masses above 42.5 M${}_{\text{Jup}}$ likely form in the same way as stars, through molecular cloud fragmentation. Troup et al. (2016) find 14 candidate BD companions around metal-poor stars ([Fe/H] $<$ -0.5), also challenging the core-accretion formation model for BDs, which they find unsurprising given that Carney et al. (2003) observed that lower metallicity populations generally have a higher number of binaries.

Here we present new BD candidate companions using data from the Sloan Digital Sky Survey III (SDSS-III) Multi-object APO Radial Velocity Exoplanet Large-area Survey (MARVELS; Ge et al. (2008)). MARVELS observed $\sim$3300 stars with $\geq$10 RV measurements. Unlike APOGEE, MARVELS primarily targeted solar-like dwarfs, and thus, is complementary to the work done by APOGEE, which primarily targets giant stars. The latest University of Florida Two Dimensional (UF2D) pipeline has RV precisions of $\sim$60-100 m s*-1* (Thomas, 2015). With this new pipeline we find 12 new low-mass companion candidates around main-sequence stars, 10 of which are in the BD mass regime.

In $\S$ 2 we discuss the detailed processing of the MARVELS data required to obtain RV measurements, including the history of its data pipelines. In $\S$ 3 we describe our companion vetting procedure designed to minimize false positive selection. In $\S$ 4 we describe our procedure for obtaining host star parameters. In $\S$ 5 we present our results and conclusions including our brown dwarf companion parameters.

2 MARVELS DATA

MARVELS was in operation from 2008 - 2012 as a part of the 6-year SDSS-III program (Eisenstein et al., 2011), a continuation of the Sloan Digital Sky Survey (SDSS; York et al., 2000). The goal of MARVELS is to detect a large sample of short-to-intermediate period giant planets obtained from a large, well-characterized, and homogeneous survey of stars with known properties. Radial velocity (RV) observations were taken using the SDSS 2.5 m telescope at Apache Point Observatory (Gunn et al., 2006). MARVELS measures RVs using a fiber-fed dispersed fixed delay interferometer (DFDI) with a moderate dispersion spectrograph. Specifically, this set-up combines a Michelson interferometer and a medium resolution (R $\sim$ 11,000) spectrograph, which overlays interferometer fringes on a long-slit stellar spectrum (Ge et al., 2009; Eisenstein et al., 2011). Unlike echelle spectrographs, which require many orders and a large detector area to achieve high RV precision, this unique setup allows simultaneous observation of 60 stars with moderate to high RV precision.

van Eyken et al. (2010) presents a complete overview of the principles and physics behind the DFDI method. The general creation process of the MARVELS DFDI spectra is explained by Ge (2002). Starlight from the telescope or light from a calibration source, a Tungsten lamp shining through a temperature-stabilized Iodine gas cell (TIO), is fed into an optical fiber. The incoming beam is split in two by a Michelson-type interferometer. A parallel glass plate creates a fixed optical delay between the two plates of the interferometer. A virtual white fringe pattern is created and re-imaged onto the entrance slit of the spectrometer, and is then dispersed by the medium-resolution grating (Ge, 2002). The dispersed fringes are recorded on a detector. The resulting fringes overlaid on the 2D spectra contain high spatial frequency Doppler information allowing for precise RV measurements. Ge et al. (2006) displayed the ability of the DFDI method with the detection of a giant planet of minimum mass $\sim$0.488 M${}_{\text{Jup}}$ (HD 102195b) using the precursor to the MARVELS instrument, the Exoplanet Tracker.

Early analysis of MARVELS data showed good instrumental performance over timescales of less than a month. By observing a RV stable star or a known planet hosting star, Ge et al. (2009) found the root mean square (RMS) errors vary between only 1-2 times the photon noise for RV measurements within one month. However, when analyzing reference stars over multiple month periods, Ge et al. (2009) found RMS errors as large as a factor of 5 times higher than the photon limiting error. They acknowledged the need for further consideration of long-term effects such as distortion, illumination, and large barycentric velocity. During MARVELS data pipeline development, Lee et al. (2011) found systematic RV errors in the form of month-to-month offsets at the level of tens to hundreds of m s*-1*. Lee et al. (2011) inferred that the systematic errors may be due to imperfections in the detailed preprocessing of the images.

Alam et al. (2015) provide an overview of the history of the MARVELS instrument and the associated data reduction. The original pipeline results are presented in the SDSS Data Release 11 (DR11) as the “cross correlation function" (CCF) and DFDI reductions, which revealed the large instrumental errors and the requirement of a new analysis method. This discovery prompted the creation of a new refurbished pipeline designed to understand and correct the systematic errors (Thomas et al., 2016). This refurbished pipeline is presented in the SDSS Data Release 12 (DR12) as the “University of Florida One Dimensional" (UF1D) reductions. Alam et al. (2015) summarize the differences between the CCF+DFDI and UF1D reductions including changes in the calibration of the wavelength solution on the detector, identification and extraction of each spectrum, and the measurement of the slant of the interferometric comb and of the resulting interference pattern of the absorption-line features. The DR12 UF1D reductions do not show great improvements over the DR11 CCF+DFDI reductions; however, the UF1D reductions resulted in a better understanding of the systematic errors of the instrument and thus allowed more accurate error estimations for individual targets.

Thomas et al. (2016) improved upon the DR12 UF1D results by adopting a ‘Lucky RV’ method, where statistically preferred regions of the spectra are used in the determination of overall RV and poorer regions are discarded once performance deteriorates for a given target (each star goes through $\sim$1,000 iterations to randomly select a specific combination and number of regions that create the best RV performance). Using this new method combined with two new filters, which result in an overall rejection of $\sim$20% of all observations, Thomas et al. (2016) was able to greatly reduce large long-term RV errors that plagued previous pipelines. The improved UF1D pipeline with Lucky RV and filters results in a RV precision of $\sim$66-114 m s*-1* (photon limits $\sim$27-90 m s*-1*) for 8-12 magnitude stable stars, respectively. The UF1D pipeline reductions identified that optical fiber illumination changes, likely isolated to fiber couplings in the gang connector, caused variations between the stellar and calibration fiber illuminations over the lifetime of the survey, causing the largest measurement error in the MARVELS RVs. However, the final UF1D pipeline still does not completely remove the systematic offsets from all targets, with an overall influence from offset errors having a median of 64.3 m s*-1*; the Lucky RV method also has the disadvantage of excluding $\sim$20% of observations Thomas et al. (2016).

After the refinement of the UF1D pipeline, a new “University of Florida Two Dimensional" (UF2D) pipeline was developed in order to use the interferometric power of MARVELS and take full advantage of the instrument’s capabilities. Rather than relying on the underlying physics of interferometry, as previous 2D DFDI pipelines generally did (Ge, 2002; van Eyken et al., 2010; Wang et al., 2011), the UF2D pipeline uses a simpler and more visually understandable method, which is outlined by Thomas (2015) and Wang, J. et al. (2012). Two new innovations, specifically 2D modeling of the underlying fringe density distortion and the correction of temporal distortions, have led to current RV performance of $\sim$60-96 m s*-1* for 8-12 magnitude RV stable stars, respectively. Figure 1 shows UF2D RV results for stable stars of various magnitudes. However, the main improvement of the UF2D pipeline is the lack of large errors found over months and years that have previously plagued the survey. Most notably, combining the results from adjacent beams of the same star improves RV accuracy, which was not the case with 1D data due to systematic errors dominating both beams. Figure 2 displays the overall RV performance for each MARVELS pipeline. We use these new UF2D RV results for our current search of substellar companions.

3 Candidate Vetting Procedure

The MARVELS survey has already identified two low-mass stellar companions to solar-type stars (MARVELS-2b, Fleming et al. (2012); MARVELS-3b, Wisniewski et al. (2012)), as well as three BD companions (MARVELS-4b, Ma et al. (2013); MARVELS-5b, Jiang et al. (2013); MARVELS-6b, De Lee et al. (2013)). Ma et al. (2016) announced the discovery of a system with both a giant planet (MARVELS-7b) and BD (MARVELS-7c) orbiting the primary star of a close binary system. The MARVELS Pilot Project also showed success in this domain with the detection of a BD (Fleming et al., 2010) and two new eclipsing binaries (Fleming et al., 2011). However, the MARVELS survey has had difficulties in this area as well. Lee et al. (2011) announced the discovery of a $\sim$28.0 M${}_{\text{Jup}}$ companion (MARVELS-1b) to TYC 1240-00945-1, but this was later determined to be a false positive in the form of a face-on double-lined binary star (Wright et al., 2013). Mack et al. (2013) warns of a similar situation where the MARVELS target TYC 3010-01494-1, a highly eccentric spectroscopic binary, initially appeared as a single star with a substellar companion with minimum mass of $\sim$50 M${}_{\text{Jup}}$ despite an extensive amount of observation and analysis.

MARVELS observed $\sim$5,500 stars during its 4-year survey, and produced 3,257 stars with 10 or more valid RV observations in our latest UF2D pipeline. Previous results show that finding substellar companions in this dataset is susceptible to error. Examining each star can also be extremely time consuming and will likely lead to biases due to expectations of the person examining the individual target. We created an automatic method to search for substellar companions that would minimize bias and time spent on each target. We focus on finding companions around solar-type stars (dwarfs/subgiants) and apply several selection criteria to avoid false positives.

3.1 Defining a Stellar Sample

We determine a giant or dwarf/subgiant designation for each star by a $J$-band reduced proper motion (RPMJ) cut, detailed in Paegert et al. (2015). Here we use the same method as Paegert et al. (2015) and compute RPMJ as follows:

[TABLE]

[TABLE]

where $J$ is the star’s 2MASS Survey (Skrutskie et al., 2006) $J$-band magnitude and $\mu_{r}$, $\mu_{d}$, and $d$ are Guide Star Catolog 2.3 (GSC 2.3; Lasker et al., 2008) proper motions in right ascension and declination (in arc seconds per year) and declination, respectively. An empirical RPMJ cut described in Collier Cameron et al. (2007) is applied:

[TABLE]

where $H$ is the star’s 2MASS Survey $H$-band magnitude.

Stars with RPMJ $\leq$ $y$ are regarded as RPMJ-dwarfs and stars with RPMJ $>$ $y$ as RPMJ-giants. With MARVELS data this method has a giant contamination rate of $\sim$4% and subgiants are mixed in with the “dwarf" sample at a level of 20%-40% (Paegert et al., 2015). Figure 3 shows the RPMJ cut on our sample of 3,257 stars with $\geq$10 RV observations, resulting in 2,378 stars classified as dwarf/subgiant. Our final stellar sample consists of 2,340 stars that have $\geq$10 RV observations, are RPMJ dwarf/subgiants, and have a time baseline longer than 300 days. These 2,340 stars have a mean baseline of $\sim$582 days and a median of 26 RV observations.

3.2 Initial Companion Selection

We first analyze our sample of 2,340 stars by determining if they are RV variable using the criteria defined in Troup et al. (2016):

[TABLE]

where v and $\sigma_{\textbf{v}}$ are the RV measurements and their uncertainties, and $\widetilde{\text{v}}$ is the median RV measurement for the star. 653 stars in our sample are RV variable with these criteria. We find a False Alarm Probability (FAP) of $\sim$0.2% for these 653 stars having real RV variation. We find a FAP by simulating RVs for each star based on a Gaussian distribution with a sigma equal to the mean RV error for that star. We calculate $\Sigma_{RV}$ for each star and find the number of simulated stars with $\Sigma_{RV}$ $\geq$ 2.5. We run this simulation 100 times and find a mean of 1.3 stars with $\Sigma_{RV}$ $\geq$ 2.5, corresponding to a FAP of (1.3/653)$\sim$0.2%. The MARVELS UF2D data has $\Sigma_{RV}$ values that peak at $\sim$0.6, which is lower than the expected value of 1 given correct error bars, suggesting a rather conservative estimate of our errors. The UF2D RV error designation for this study is similar to the UF1D pipeline error designation, detailed in (Thomas et al., 2016). The RV errors are based on photon noise, general pipeline errors, and illumination profile (or Y-profile, see Thomas et al. (2016)) variation errors. We find that these errors are somewhat conservative; however, using more conservative errors ensures that only stars with strong RV signals will be considered as new candidates for this study.

We search for periodic signals in these 653 RV variable stars using a Lomb-Scargle (LS) periodogram (Lomb, 1976; Scargle, 1982), which is well suited to find periodic signals in unevenly spaced RV data. For each RV variable star we find the period with the max LS periodogram power between one and 750 days ($\sim$max temporal baseline of MARVELS data), which we use as an initial guess for our Keplerian fitting. We also fit a basic linear trend to each star’s RV data, subtract this trend, and find the max LS periodogram power for this data. For Keplerian fitting we use the RV orbit fitting software, RVLIN, an IDL111IDL is a commercial programming language and environment by Harris Geospatial Solutions: http://www.harrisgeospatial.com/IDL package developed by Wright & Howard (2012). This software uses a Levenberg-Marquardt algorithm to fit RV data and find the best-fitting orbital parameters. Each RV variable star has eight initial fits. We input four different eccentricities (0.05, 0.3, 0.6, and 0.9) for both of the max LS periodogram periods found with and without a linear trend removed. We include a linear component in the Keplerian fit for the max LS periodogram periods found after removing a linear trend in the RV.

All fits are run through numerous selection criteria designed to eliminate poor fits to the data. We use several criteria from Troup et al. (2016), who find 382 stellar and substellar companions in the APOGEE data. Some of these criteria use $U_{N}$, the uniformity index (Madore & Freedman, 2005), which is used to quantify the phase coverage of the data; as well as $V_{N}$, the “velocity" uniformity index (Troup et al., 2016), which quantifies the velocity coverage and prevents a selection program from settling on a very eccentric orbit that is unsupported by the data. $U_{N}$ is defined as:

[TABLE]

where the values $\phi_{i}$ are the sorted phases associated with the corresponding Modified Julian Date (MJD) of the measurement i, and $\phi_{N+1}$ = $\phi_{1}$ + 1. $U_{N}$ is normalized such that $0\leq U_{N}\leq 1$, where $U_{N}$ = 1 indicates a curve evenly sampled in phase space. $V_{N}$ is defined as:

[TABLE]

where $\nu_{i}$ = (v${}_{\text{i}}$ - v${}_{\text{min}}$)/(v${}_{\text{max}}$ - v${}_{\text{min}}$). The values of v${}_{\text{i}}$ are the radial velocity measurements, sorted by their value. If a fit has a linear trend applied, then that trend is subtracted from the v${}_{\text{i}}$ values. For the fits with no linear trend applied, we subtract the average of the raw velocities: v${}_{\text{i}}$ = v${}_{\text{raw,i}}$ - $\bar{\text{v}}_{\text{raw}}$. We also make use of the modified $\chi^{2}$ statistic, which is defined by Troup et al. (2016):

[TABLE]

where $\chi^{2}_{\text{red}}$ is the traditional reduced $\chi^{2}$ goodness-of-fit statistic.

A fit must meet these 11 criteria to be considered valid:

1. (i)

Periods cannot be within 5% of 3, 2, 1, 1/2, or 1/3

type– days (to avoid one-day aliased systems) 2. (ii)

Eccentricities cannot be extremely high ($e>0.934$) 3. (iii)

Fits must have good phase and velocity coverage:

type— ($U_{N}V_{N}\geq 0.5$) 4. (iv)

$\frac{K}{\widetilde{\sigma}_{\text{v}}}\geq 3+3(1-V_{N})e$ 5. (v)

$\frac{K}{|\Delta\widetilde{\text{v}}|}\geq 3+3(1-V_{N})e$ 6. (vi)

$\chi^{2}_{\text{mod}}\leq\frac{K/\widetilde{\sigma}_{\text{v}}}{3+3(1-V_{N})e}$ 7. (vii)

$\chi^{2}_{\text{mod}}\leq\frac{K/|\Delta\widetilde{\text{v}}|}{3+3(1-V_{N})e}$ 8. (viii)

$\frac{K}{\mbox{stddev}(|\textbf{v}-\widetilde{\text{v}}|)}<4$ 9. (ix)

If $e$ $>$ 0.4 then there must be an RV point in the

type— phase space ranges $\phi$ $<$ 0.1 and $\phi$ $>$ 0.9. 10. (x)

$\Delta T\geq 2P$ 11. (xi)

$\frac{K}{\text{max}(\textbf{v})-\text{min}(\textbf{v})}\geq 1/3$

where $\widetilde{\text{v}}$ is the median RV measurement for the star, $\widetilde{\sigma}_{\text{v}}$ is the median RV uncertainty, $K$ is the RV semi-amplitude of the best-fit model, $|\Delta\widetilde{\text{v}}|$ is the median absolute residuals of the best-fit model, v are the RV measurements, $\Delta T$ is the temporal baseline of the star, P the period of the fit, $U_{N}$ is the uniformity index, $V_{N}$ is the velocity uniformity index, and $\chi^{2}_{\text{mod}}$ is the modified $\chi^{2}$ statistic.

Criteria 1 - 7 are detailed in Troup et al. (2016). Criteria 8 - 11 are additional criteria to those in Troup et al. (2016), which further eliminate poor fits to the RV data. Criteria 8 and 9 are designed to eliminate highly eccentric fits with large RV semi-amplitudes that are clearly not supported by the data. We apply a strict period coverage with criterion 10, which ensures that the RV data covers at least two full phases of the fitted period. We eliminate unlikely fits with linear trends much larger than the semi-amplitude using criterion 11.

Of the initial 653 RV variable MARVELS stars, 84 have at least one Keplerian fit that meets all 11 selection criteria. If a star has multiple Keplerian fits that match all 11 criteria, we use the fit with the lowest $\chi^{2}_{\text{mod}}$ and no linear trend component as our fit for the remainder of the vetting procedure. We use a Keplerian fit with a linear trend component if only fits with a linear trend passed all the criteria. Assuming a 1 M*⊙* host star, 27 of these companion candidates have M sini values within 10 M${}_{\text{Jup}}$ - 110 M${}_{\text{Jup}}$. We use these 27 candidates as our final BD candidates from this initial selection process. In $\S$ A we determine rough dependencies on our results for each of the 11 criteria. Given the low BD occurrence rates of previous studies of $\lesssim$1% (Vogt et al., 2002; Patel et al., 2007; Wittenmyer et al., 2009; Sahlmann et al., 2011; Santerne et al., 2016), we would expect to find $\lesssim$ 23 BDs in our stellar sample. Occurrence rate expectations combined with prior studies, which exhibit the likelihood of false positives in the MARVELS data, deems further analysis necessary to validate these candidates.

3.3 Periodogram False Alarm Probability

A Lomb-Scargle (LS) periodogram is well suited to find periodic signals in unevenly spaced RV data, but periodogram peaks do not necessarily signify real substellar companions. Although the periodogram contains information on how well an oscillating signal at any grid frequency fits the observed data, it also includes signs of other undesirable effects including aliases, parasite frequencies, and spectral leakage (Suveges, 2012). The false alarm probability (FAP) evaluates the statistical significance of a periodogram peak by testing the zero hypothesis of the observed time series (pure white noise) against a periodic deterministic signal.

We assess the FAP for a given orbit fit by using a randomized bootstrapping method to compare how the likelihood of the observed RVs compare to a random set with the same sampling (Suveges, 2012). We shuffle all of the RV results among the existing observation times and then record the peak power in the LS periodogram from this scrambled data set. This is repeated for 10,000 iterations. If more than 0.1% of the peak powers in the randomized data exceed the power of the model period fit (0.1% FAP), then we consider the candidate to be a likely false positive. However, if a companion’s model fit has an eccentricity $>$ 0.5, the model period’s peak must only be within 90% of the max peak of the bootstrapped data. If a candidate’s fit has a linear trend, we first remove the linear trend from the RV before computing the LS periodogram values. Figures 11 - 22 show periodograms for each of the final BD candidates. This test eliminates four of the potential BD candidates, leaving a total of 23 BD candidates. A less strenuous 1% FAP requirement still excludes these four candidates.

Finally we analyze these remaining 23 candidates by looking at both their RV measurements and LS periodograms. We find that our period fit for TYC 1717-02014-1 sits just above the 0.1% FAP level, and when looking at the RV, a stellar companion with a long period seems just as likely as the low-mass fit our vetting procedure produces (figure 4). Thus, we remove TYC 1717-02014-1 from our potential BD candidates. We also find that three of our potential BD candidates all have similar RV fits and LS periodogram peaks (figure 5). These three candidates are all on the same SDSS-III MARVELS plate (each plate contains fiber holes for the 60 stars that were observed at the same time) and may suffer from the same systematic errors associated with the measurements of this plate. We remove TYC 3545-02661-1, TYC 3549-00411-1, and TYC 3545-00813-1 from our list of potential BD companions as these RV variations may be due to systematic errors rather than real low-mass companions. This final periodogram analysis leaves us with 19 potential BD candidates.

3.4 Stellar Contamination

RV signals containing amplitudes compatible with substellar companions are not necessarily planets or BDs; the RV signal might be mimicked by stellar activity (e.g. Queloz et al., 2001). Wright et al. (2013) and Mack et al. (2013) displayed this type of contamination in previous MARVELS companion candidates. Many methods to reduce this type of contamination employ the analysis of a cross-correlation function (CCF), which can find a star’s RV by measuring the Doppler shift of spectral lines. The CCF technique averages hundreds to thousands of spectral lines to create a spectral line profile that can be approximated by a Gaussian function. The CCF has a much higher SNR compared with individual lines and can be characterized by a given location in the RV space, a contrast, and a FWHM. An asymmetry in the observed line profile can be produced by the convective blue shift, pulsations, stellar spots, and contaminating stars (Santerne et al., 2015).

We use the bisector velocity span (BIS; Queloz et al., 2001) to measure line-profile asymmetry in the CCFs of our potential BD candidates. The BIS consists of computing the velocity span between the average of the top portion of the bisector and the bottom portion of the bisector, where the bisector is defined as the mid-points of horizontal line segments bounded by the sides of the line profile. We use the definition and limits described in Santerne et al. (2015) for our BIS calculations. We also use the FWHM of the CCF line profiles to measure CCF line-profile asymmetry. The FWHM may present a clear correlation with RV when other line-profile asymmetry diagnoses do not; therefore, the FWHM is a good diagnosis to use as a complement to other line-asymmetry diagnoses (Santerne et al., 2015). For each BD companion candidate we calculated the Kendall rank correlation coefficient between their RV and BIS and RV and FWHM values. However, as pointed out by Mack et al. (2013), MARVELS low resolution makes line bisector analysis less reliable than other studies which use high resolution spectra for this method. Thus, we analyzed both the CCF line profiles and kendall rank correlation p-values for these candidates. For our CCF line profiles we take the average CCF line profile of both spectral beams of the star for each observation. The CCF line profile may contain an obvious second component, as seen in figure 6a, which shows the CCF line profiles and line-asymmetry diagnoses plots for the double-lined spectroscopic binary TYC 3010-01494-1 (Mack et al., 2013).

The 19 BD candidates from the section above includes TYC 1240-00945-1, the face-on binary previously reported to be a BD (Lee et al., 2011; Wright et al., 2013). Figure 6b shows the RV and line-asymmetry diagnoses plots with the associated p-values and the CCF line profiles for this star. We see that there is not a visible second component in the CCF line profile but this star has a low BIS p-value of $\sim$0.003. In order to set a significance level, we calculate p-values for the $\sim$2,340 MARVELS stars with $\geq$ 10 RV observations, time baselines longer than 300 days, and RPMJ dwarf/subgiant designation. This sample has a mean RV and BIS correlation p-value of $\sim$0.31 and a mean RV and FWHM correlation p-value of $\sim$0.30. The BIS p-value of 0.00265 for the known binary TYC 1240-00945-1 is greater than $\sim$16% of stars in this sample, while the relatively high FWHM p-value of $\sim$0.54 for this star is greater than $\sim$75% of stars in this sample. We set our p-value significance levels at the bottom 20% of this sample, which corresponds to a BIS p-value of $\sim$0.008 or lower and a FWHM p-value of $\sim$0.01 or lower.

One other star besides TYC 1240-00945-1 has a BIS p-value below 0.008: TYC 3005-00460-1. Considering the CCF line profile variation and low p-value seen in figure 7, as well as no associated linear trend in the RV fit that may account for such a variation, we remove TYC 3005-00460-1 as a BD candidate. Two stars have a FWHM p-value below 0.01: TYC 3413-02367-1 (MCKGS1-94 from Ghezzi et al. (2014)) and TYC 5634-00322-1. We also remove these stars as BD candidates. Figures 11 - 22 show the CCF line profiles and line-asymmetry diagnoses plots for all of the new MARVELS BD candidates. TYC 0173-02410-1 (MARVELS-9b, figure 12) shows relatively large CCF line profile variation, but our RV fit contains a linear component that may account for the CCF asymmetry.

3.5 Previous Candidates & Vetting Results

We removed TYC 1240-945-1 (MARVELS-1b, Lee et al., 2011; Wright et al., 2013) as well as TYC 3005-00460-1, TYC 3413-02367-1 (MCKGS1-94), and TYC 5634-00322-1 in the previous section, leaving 15 BD candidates left in our vetting procedure. However, three of these 15 candidates were previously published: the $\sim$32 M${}_{\text{Jup}}$ BD around GSC 03546-01452 (MARVELS-6b, De Lee et al., 2013), the $\sim$65 M${}_{\text{Jup}}$ BD around HIP 67526 (MARVELS-5b, Jiang et al., 2013), and the $\sim$98 M${}_{\text{Jup}}$ stellar companion around TYC 4110-01037-1 (MARVELS-3b, Wisniewski et al., 2012). We find similar companion parameters as previous studies for these three stars using the new UF2D data. With the exclusion of these three companions, our full vetting process produces 12 new BD candidates.

The other MARVELS stars with published companions or double-lined spectroscopic binaries (SB2s) did not make it through the vetting procedure for various reasons. We did not find a companion around TYC 2087-00255-1 (MARVELS-4b, Ma et al., 2013) due to three outliers in the new UF2D data, which prevent the vetting procedure to converge on the best fit initially. After removing these three outliers, we recover the companion reported in Ma et al. (2013). We did not find any companions around TYC 2930-00872-1 (MARVELS-2b, Fleming et al., 2012); the initial fit was poor as the star likely requires a multiple companion fit. The vetting criteria also eliminated HD 87646 (MARVELS-7b,c, Ma et al., 2016), again likely due to the multiple companion fitting requirements. The UF2D pipeline produces large outliers for the SB2 TYC 3010-01494-1 (Mack et al., 2013), preventing a good-fit and possible false-positive that could have made it through all of our initial criteria. The SB2 TYC 1240-00945-1 (MARVELS-1b, Lee et al., 2011; Wright et al., 2013) has a companion fit that makes it through all 11 initial criteria, but is excluded in the stellar contamination section (see figure 6b).

4 Host Star Parameters

Ghezzi et al. (2014) determine precise atmospheric parameters (effective temperature T${}_{\text{eff}}$, metallicity [Fe/H], and surface gravity log g) using a method based on spectral indices. This method is necessary because blending of the atomic lines is severe at the moderate resolution of the MARVELS spectra, precluding individual line analysis through classical spectroscopic methods (e.g., excitation and ionization equilibria). The spectral indices method does not rely on the accuracy of single atomic lines that can be blended together in moderate and low resolution spectra. Ghezzi et al. (2014) test the spectral indices method using 30 MARVELS stars that also have reliable atmospheric parameters derived from high-resolution spectra and spectroscopic analysis based on the excitation and ionization equilibria method. Using MARVELS spectra and the spectral indices method, Ghezzi et al. (2014) recover the parameters within 80 K for T${}_{\text{eff}}$, 0.05 dex for [Fe/H], and 0.15 dex for log g. In a separate test using 138 stars from the ELODIE stellar library, they recover the literature atmospheric parameters within 125 K for T${}_{\text{eff}}$, 0.10 dex for [Fe/H], and 0.29 dex for log g.

Of the 12 new companions, four were previously considered to be MARVELS candidates, for which Ghezzi et al. (2014) obtained high-resolution spectra and analysis. We use Ghezzi et al. (2014)’s high resolution values for these four stars, while the other stellar parameters are found with the spectral indices method using the same pipeline developed by Ghezzi et al. (2014). However, we use input spectra that were pre-processed using the newest UF1D/2D pipeline (Thomas et al., 2016; Thomas, 2015), while Ghezzi et al. (2014) used input spectra that were pre-processed using the older CCF+DFDI MARVELS pipeline released in the SDSS DR11. We tested the same 30 MARVELS stars as Ghezzi et al. (2014) using the newest UF1D input spectra and recovered parameters within 85 K for T${}_{\text{eff}}$, 0.06 dex for [Fe/H], and 0.17 dex for log g compared to the high-resolution analysis results. Atmospheric parameters and their associated uncertainties for the indices method are found with a reduced $\chi^{2}$ minimization process, which is detailed in Section 4.4 of Ghezzi et al. (2014). The T${}_{\text{eff}}$, [Fe/H], and log g values for the companion host stars are presented in Table 1.

We determine the mass and radius of each host star from T${}_{\text{eff}}$, [Fe/H], and log g using the empirical polynomial relations of Torres et al. (2010), which were derived from a sample of eclipsing binaries with precisely measured masses and radii. We estimate the uncertainties in M*⋆* and R*⋆* by propagating the uncertainties in T${}_{\text{eff}}$, [Fe/H], and log g using the covariance matrices of the Torres et al. (2010) relations (kindly provided by G. Torres). The mass and radius of each host star is also presented in Table 1.

5 RESULTS & Conclusions

5.1 Brown Dwarf Parameters

Our candidate filtering process produced 12 new companions within $\sim$10 - 110 M${}_{\text{Jup}}$. We run these stars through RVLIN a final time to obtain an optimal fit. Figures 11 - 22 show the best fit to the RV data for each of the new MARVELS companion candidates. The orbital parameters for our final fits of the MARVELS BD candidates introduced in this paper are presented in Table 2. Values for M sini are determined using the host star mass derived in the previous section. Uncertainties in all parameters are calculated through bootstrapping (with 1,000 bootstrap replicates) using the publicly available BOOTTRAN package, which is detailed in Wang, X. et al. (2012). Our final parameters produce two low-mass stellar companions (M sini $\sim$97.6 and 100.2 M${}_{\text{Jup}}$), and 10 new BD candidates in the mass range 13 M${}_{\text{Jup}}$ $<$ M sini $<$ 80 M${}_{\text{Jup}}$. TYC 3547-01007-1 (MARVELS-17b, figure 20) resides just above the giant planet regime with M sini = 13.1 $\pm$ 1.6 M${}_{\text{Jup}}$. Given the errors of M sini on MARVELS-17b it may be a giant planet if its orbit is edge-on.

5.2 Caveats

We designed our vetting procedure to only extract substellar companion candidates with a low likelihood of being false positives. However, the new substellar companions presented in this work are still only candidates. Follow-up observations are required to confirm them as true BDs. Some follow-up observations have already begun using the EXtremely high Precision Extrasolar Planet Tracker (EXPERT; Ge et al., 2010) instrument at the Kitt Peak National Observatory. Additional CCF line profile analysis with high resolution spectra would help confirm that these candidates are not double-lined spectroscopic binaries (SB2s).

Some candidates pass all of our selection criteria but have caveats which make them more prone to being false positives. These caveats consist of weak phase coverage, multiple periodogram peaks, and aliases that make other RV fits just as likely as those presented. More observations are required to fully constrain some of our new candidate parameters.

Specifically we note the poor phase coverage of MARVELS-12b (figure 15). Given its high eccentricity, this candidate could have failed to pass criterion 9; however, this star had at least one fit pass all 11 initial criteria. When fitting some stars a final time with RVLIN we may find a better fit with higher eccentricity, which also passes the periodogram FAP test. More phase coverage is needed to validate this candidate, but given that at least one fit can pass all 11 criteria and it passes both our stellar contamination and periodogram tests we put a high likelihood on this star having some type of low-mass companion.

5.3 Statistical Results & Conclusions

Ma & Ge (2014) collected data from literature about currently known BD candidate companions around FGK-type stars. While examining the period-mass distribution of BDs they found a statistically significant gap in the short period (P $<$ 100 d) and medium-mass region (35 M${}_{\text{Jup}}$ $<$ M${}_{\text{BD}}$ $<$ 55 M${}_{\text{Jup}}$), which they deemed the ‘driest’ part of the brown dwarf desert and likely to be real as it is unlikely caused by detection sensitivity (RV precision) or observation biases (due to survey incompleteness). Three of our new BD candidates, however, are within this driest part of the desert. Figure 8 shows the M sini vs period plot of Ma & Ge (2014), but now with candidates from Table 7 of Wilson et al. (2016) combined with our new MARVELS candidates.

As previously discussed, Ma & Ge (2014) concluded that the orbital eccentricities for BDs with masses lower than 42.5 M${}_{\text{Jup}}$ and those higher than 42.5 M${}_{\text{Jup}}$ have noticeable differences in their orbital eccentricities. Their study showed that BDs below $\sim$42.5 M${}_{\text{Jup}}$ trend such that more massive BDs tend to have lower maximum eccentricities, similar to giant planets, while BDs above this mass threshold do not show such a trend, instead showing more diversity in their eccentricities. This trend is consistent with less massive BD formation occurring in protoplanetary discs and subsequential scattering with other objects in the disc pumping these BDs to higher eccentricities. In the scattering model (Chatterjee et al., 2008; Ford & Rasio, 2008), it is harder to pump a higher mass BD to higher eccentricity. Figure 9 shows the eccentricity vs M sini plot from Ma & Ge (2014), now with Wilson et al. (2016) BDs and new MARVELS BD candidates as in figure 8. The two violet curves show eccentricities if BDs were scattered by another object of 20 or 25 M${}_{\text{Jup}}$. We see a few outliers from Wilson et al. (2016) as well as one new MARVELS BD candidate in the high eccentricity region above the Ford & Rasio (2008) curves; however, in general our new BD candidates still support a two population trend.

Two of our candidates have host stars that are classified as subgiants according to Troup et al. (2016) criteria: T${}_{\text{eff}}$ $>$ 4800 K, log g $\geq$ 3.7 + 0.1[Fe/H] if T${}_{\text{eff}}$ < 5500 K, and log g $\leq$ 4 - (7 x 10*-5*)(T${}_{\text{eff}}$ - 8000 K); the third criterion was derived by mapping the main-sequence turnoff in Dartmoth isochrones. This includes TYC 3097-00398-1 (MARVELS-15, log g $\sim$ 3.85 dex) which hosts the $\sim$98 M${}_{\text{Jup}}$ stellar companion (figure 18) and TYC 3469-00492-1 (MARVELS-12, log g $\sim$ 3.78 dex) which hosts the $\sim$33 M${}_{\text{Jup}}$ BD candidate with a $\sim$253 day period (figure 15).

The host stars of our 10 new BD candidates (M sini $\leq$ 80 M${}_{\text{Jup}}$) have a mean metallicity of [Fe/H] = 0.04 $\pm$ 0.07 dex. When including MARVELS-4b, MARVELS-5b, MARVELS-6b, and MARVELS-7c with our new candidates we find a mean metallicity of [Fe/H] = 0.03 $\pm$ 0.08 dex. This result agrees with Mata Sánchez et al. (2014) who study chemical abundances of stars with BD companions. Mata Sánchez et al. (2014) find that stars with BD companions have $\alpha$-elements and iron peak elements at about solar abundance, whereas stars with low-mass and high-mass planets have [$X_{\alpha}$/H] and [$X_{\text{Fe}}$/H] peak abundaces of $\sim$0.1 dex and $\sim$0.15 dex, respectively. Three of the 10 new BD candidates have host stars with subsolar metallicity. Core accretion, the giant planet formation mechanism dependent on high metallicity (Ida & Lin, 2004), is less likely than other mechanisms for these BD candidates.

In this work we sifted through MARVELS data to create a sample of $\sim$2,340 dwarf/subgiant stars with reasonable RV coverage ($\geq$10 observations) and time baselines ($\geq$ 300 days) and found 10 new BD candidates. Thus, including the three previously found MARVELS BDs with periods less than 300 days (MARVELS-4b, Ma et al. (2013); MARVELS-5b, Jiang et al. (2013); MARVELS-6b, De Lee et al. (2013)), we estimate the BD occurrence rate for solar-type stars up to $\sim$300 days (approximately half the typical baseline of our stellar sample) to be $\sim$0.56%. Unlike Troup et al. (2016), our results do not show a truncated brown dwarf desert, but one that extends up to at least 1 AU. However, Troup et al. (2016)’s APOGEE sample is dominated by evolved stars, whereas our stellar sample was designed to exclude giant stars. Our work gives further support for the brown dwarf desert around solar-type stars.

Acknowledgements

Funding for the MARVELS multi-object Doppler instrument was provided by the W. M. Keck Foundation and NSF with grant AST-0705139. The MARVELS survey was partially funded by the SDSS-III consortium, NSF grant AST-0705139, NASA with grant NNX07AP14G, and the University of Florida. Funding for SDSS-III has been provided by the Alfred P. Sloan Foundation, the Participating Institutions, the National Science Foundation, and the U.S. Department of Energy Office of Science. L.G. would like to thank the financial support provided by CAPES/PNPD Fellowship. G.F.P.M. acknowledges the CNPq research grant 474972/2009-7.

Appendix A Vetting Criteria Dependencies

Here we determine rough estimates for the dependence of our results on each vetting criterion defined in $\S$ 3.2. For each criterion we add and subtract 1$\sigma$ from the main dependence and determine the number of stars that have companion fits that pass all 11 criteria, labeled as # pass in table 3, without changing any of the other 10 criteria. Out of these stars with passing companion fits we then find the number of companions that are within the mass range we consider for BD candidates (M sini values within 10 M${}_{\text{Jup}}$ - 110 M${}_{\text{Jup}}$ assuming a 1 M*⊙* host star), labeled as #BD pass in table 3. The original criteria have 84 stars with companion fits that pass all 11 criteria and 27 of them are within our mass range for BD candidate consideration. Figure 10 further illustrates criteria dependencies and shows the results of each BD candidate as we vary each criterion -1$\sigma$ and +1$\sigma$. As illustrated in both figure 10 and table 3, the final number of candidates is most sensitive to the reduction of criterion 8. This is expected because making this criterion more strict will remove several candidates with relatively high semi-amplitudes compared to their velocity variations.

Appendix B New MARVELS Substellar Companions

Here we include both the full-length (top left) and phased (top right) radial velocity plots, Lomb-Scargle periodograms (bottom left), CCF line profiles (bottom middle), and BIS and FWHM vs RV plots (bottom right) for each new MARVELS companion candidate. Radial velocity measurements are obtained from MARVELS and processed with the UF2D pipeline. The top panel for each radial velocity plot shows the best fit to the data, and the bottom panel shows the best-fit residuals. Lomb-Scargle periodograms show False Alarm Probabilities (FAPs) that are determined via a randomized bootstrapping method with 10,000 iterations. The dotted lines represent 50% FAP, dashed lines represent 1% FAP, and the dashed and dotted lines represent 0.1% FAP. Filled circles represent the best-fit period of each companion. The best-fit period must be above the 0.1% FAP, unless companion’s fit has an eccentricity $>$ 0.5, then the period’s peak must only be within 90% of the max peak of the bootstrapped data. MARVELS-9b (figure 12), MARVELS-10b (figure 13), and MARVELS-17b (figure 20) had the best-fit linear trend removed from the RV data before running the LS periodograms. The BIS and FWHM vs RV plots show the bisector velocity span (BIS) and the FWHM of the CCF line profiles as a function of radial velocity for each of the MARVELS candidate host stars. Both BIS and FWHM values have the mean value of the star subtracted before plotting. We show the p-value of the Kendall rank correlation between the BIS and RV and FWHM and RV values in each plot.