MOA-bin-29b : A Microlensing Gas Giant Planet Orbiting a Low-mass Host Star

TL;DR

This paper reports the discovery of a gas giant planet orbiting a low-mass host star in the Galactic bulge through microlensing, with detailed analysis of the system's physical properties and degeneracies.

Contribution

First detection of a gas giant orbiting a low-mass host star via microlensing, with Bayesian analysis estimating the system's physical parameters.

Findings

Gas giant mass: approximately 0.63 Jupiter masses

Host star likely a brown dwarf or late M-dwarf

System located about 6.89 kpc in the Galactic bulge

Abstract

We report the discovery of a gas giant planet orbiting a low-mass host star in the microlensing event MOA-bin-29 that occurred in 2006. We find five degenerate solutions with the planet/host-star mass ratio of $q \sim 10^{-2}$. The Einstein radius crossing time of all models are relatively short ($\sim 4-7$ days), which indicates that the mass of host star is likely low. The measured lens-source proper motion is $5-9$ ${\rm mas}\ {\rm yr}^{-1}$ depending on the models. Since only finite source effects are detected, we conduct a Bayesian analysis in order to obtain the posterior probability distribution of the lens physical properties. As a result, we find the lens system is likely to be a gas giant orbiting a brown dwarf or a very late M-dwarf in the Galactic bulge. The probability distributions of the physical parameters for the five degenerate models are consistent within the range of error. By combining these probability distributions, we conclude that the lens system is a gas giant with a mass of $M_{\rm p} = 0.63^{+1.13}_{-0.39}\ M_{\rm Jup}$ orbiting a brown dwarf with a mass of $M_{\rm h} = 0.06^{+0.11}_{-0.04}\ M_\odot$ at a projected star-planet separation of $r_\perp = 0.53^{+0.89}_{-0.18}\ {\rm au}$. The lens distance is $D_{\rm L} = 6.89^{+1.19}_{-1.19}\ {\rm kpc}$, i.e., likely within the Galactic bulge.

Tables

Table 1. Table 1: The number of data points in the light curves and the normalization parameters

Telescope	Filter	$k$	$e_{\min }$	Number of Data
MOA	$R_{\mathrm{MOA}}$	1.105	0.003	29094
OGLE	$I$	1.058	0.003	871

Table 2. Table 2: The Best-fit parameters for five competing models

Parameters	Unit	wide1	wide2	wide3	close1	close2
$t_0$	HJD-2450000	3931.097	3931.098	3931.080	3931.135	3931.139
		0.003	0.026	0.009	0.002	0.004
$t_E$	days	3.887	3.139	7.014	5.183	5.191
		0.210	0.150	0.920	0.860	0.950
$u_0$	$\times {10}^{-2}$	1.090	4.849	0.893	1.249	1.351
		0.076	0.690	0.290	0.320	0.320
$q$	$\times {10}^{-2}$	1.631	2.038	0.614	1.170	0.827
		0.240	0.630	0.210	0.240	0.210
$s$		1.164	1.219	1.745	0.537	0.589
		0.009	0.030	0.097	0.024	0.032
$\alpha$	radian	2.695	3.369	3.103	3.852	3.387
		0.027	0.073	0.048	0.04	0.076
$\rho$	$\times {10}^{-2}$	0.803	1.233	0.271	0.691	0.420
		0.092	0.120	0.097	a	a
$\Delta {\chi }^2$			5.55	7.08	13.24	13.89

Table 3. Table 3: The Best-fit parameters and the median value with 68.3% confidence interval derived from MCMC chains for the 2S1L model

Parameters	Unit	Best fit	MCMC
$t_{0,1}$	HJD-2450000	3931.151	${3931.150}_{-0.001}^{+0.001}$
$t_{0,2}$	HJD-2450000	3931.090	${3931.082}_{-0.024}^{+0.026}$
$t_E$	days	4.486	${3.998}_{-0.689}^{+0.693}$
$u_{0,1}$	$\times {10}^{-2}$	-0.652	$-{0.905}_{-0.212}^{+0.218}$
$u_{0,2}$	$\times {10}^{-2}$	4.732	${5.754}_{-1.617}^{+1.668}$
${\rho }_1$	$\times {10}^{-2}$	0.287	${0.636}_{-0.426}^{+0.442}$
${\rho }_2$	$\times {10}^{-2}$	4.788	${2.606}_{-1.552}^{+1.546}$
$q_{F,\mathrm{OGLE}}$		1.309	${12.351}_{-10.566}^{+13.266}$
$q_{F,\mathrm{MOA}}$		1.004	${0.911}_{-0.175}^{+0.165}$

Table 4. Table 4: The intrinsic Source Color, intrinsic Magnitude, Angular source Radius, Angular Einstein Radius, and the lens-source relative proper motion

Parameters	Unit	wide1	wide2	wide3	close1	close2
$I_s$	mag	19.893	19.239	20.428	20.091	20.033
		0.096	0.112	0.288	0.269	0.253
${\left(V-I\right)}_s$	mag	1.023	0.842	1.270	1.120	1.099
		0.176	0.112	0.261	0.210	0.209
${\theta }_{\ast }$	$\mu$as	0.448	0.507	0.444	0.449	0.451
		0.079	0.066	0.127	0.107	0.105
${\theta }_E$	mas	0.056	0.041	0.164	$>0.042$	$>0.068$
		0.012	0.007	0.075
${\mu }_{\mathrm{rel}}$	mas ${\mathrm{yr}}^{-1}$	5.242	4.793	8.552	$>2.807$	$>4.436$
		1.144	0.823	4.077

Table 5. Table 5: The lens physical parameters as a result of the Bayesian analysis

Parameters	Unit	wide1	wide2	wide3	close1	close2
$D_L$	kpc	${7.12}_{-1.05}^{+1.15}$	${7.23}_{-1.01}^{+1.15}$	${6.57}_{-1.31}^{+1.23}$	${7.07}_{-1.07}^{+1.15}$	${6.73}_{-1.21}^{+1.21}$
$M_L$	$M_{☉}$	${0.03}_{-0.02}^{+0.07}$	${0.02}_{-0.01}^{+0.06}$	${0.10}_{-0.06}^{+0.16}$	${0.05}_{-0.03}^{+0.08}$	${0.07}_{-0.04}^{+0.11}$
$M_p$	$M_{\mathrm{Jup}}$	${0.60}_{-0.38}^{+1.20}$	${0.48}_{-0.31}^{+1.22}$	${0.63}_{-0.36}^{+1.05}$	${0.58}_{-0.37}^{+1.04}$	${0.58}_{-0.35}^{+0.97}$
$r_{⟂}$	au	${0.48}_{-0.10}^{+0.12}$	${0.37}_{-0.07}^{+0.08}$	${1.54}_{-0.45}^{+0.52}$	${0.27}_{-0.07}^{+0.08}$	${0.39}_{-0.11}^{+0.12}$
$a_{3D}$	au	${0.59}_{-0.16}^{+0.35}$	${0.46}_{-0.11}^{+0.27}$	${1.92}_{-0.64}^{+1.18}$	${0.34}_{-0.10}^{+0.21}$	${0.49}_{-0.16}^{+0.30}$
$I$	mag	${36.35}_{-8.62}^{+3.66}$	${37.24}_{-8.57}^{+4.94}$	${27.77}_{-2.68}^{+7.77}$	${35.56}_{-8.54}^{+2.96}$	${33.59}_{-7.56}^{+3.34}$
$V$	mag	${44.61}_{-11.92}^{+2.82}$	${45.34}_{-11.60}^{+4.35}$	${32.85}_{-2.95}^{+11.00}$	${43.94}_{-12.27}^{+2.23}$	${41.62}_{-11.35}^{+3.35}$
$H$	mag	${32.34}_{-8.10}^{+4.03}$	${33.83}_{-8.99}^{+4.98}$	${24.14}_{-1.99}^{+6.99}$	${30.98}_{-7.21}^{+3.95}$	${28.58}_{-5.58}^{+4.58}$
$K$	mag	${33.11}_{-9.40}^{+7.63}$	${35.33}_{-11.14}^{+8.14}$	${23.66}_{-1.92}^{+8.07}$	${31.09}_{-7.82}^{+5.73}$	${28.29}_{-5.73}^{+6.75}$

Equations

$${\sigma }_i^{\prime }=k\sqrt{{\sigma }_i^2+e_{\min }^2},$$

$$F(t)=A(t)F_S+F_b,$$

$$\Delta F(t)$$

$$F_b^{\prime }$$

$$S_{\lambda }(\vartheta )=S_{\lambda }(0)[1-u_{\lambda }(1-\cos (\vartheta ))],$$

$$A=\frac{A_1{F}_1+A_2{F}_2}{F_1+F_2}=\frac{A_1+q_F{A}_2}{1+q_F},$$

$$R_{\mathrm{MOA}}-I_{\mathrm{OGLE}}=(-28.32\pm 0.01)+(0.16\pm 0.01)(V-I{)}_{\mathrm{OGLE}}.$$

$$(V-I,I{)}_{\mathrm{S},0}=(1.02,19.89)\pm (0.18,0.10).$$

$$\log ({\theta }_{\mathrm{LD}})=0.501414+0.419685(V-I)-0.2I,$$

Full text

gravitational lensing: micro — planets and satellites: detection

1 Introduction

More than 3900 exoplanets have been discovered since the first discovery of an exoplanet orbiting a main-sequence star in 1995 (Mayor & Queloz, 1995), including various planetary systems, such as hot Jupiters and super Jupiters. Most known exoplanets have been found by the radial velocity (Butler et al., 2006; Pepe et al., 2011) and transit methods (Batalha et al., 2013), which are most sensitive to massive planets in close orbits. Direct imaging has found young giant planets in very wide orbits.

The gravitational microlensing method has a unique planet sensitivity to planets down to low masses (Bennett & Rhie, 1996) in wide orbits, just beyond the snow line (Gould & Loeb, 1992). Exoplanet searches by using the microlensing were first proposed by Mao & Paczyński (1991), and over 90 planets have been discovered by this method to date. Gravitational microlensing occurs when a foreground lens star crosses the line of sight between an observer and a background source star by chance. The gravity of the lens star bends the light from the source star and magnifies its brightness. If the lens star has a companion, its gravity affects the magnification of the source star. The microlensing method does not depend on the brightness of the lens objects. So we can discover low-mass companions around faint and/or distant host, such as M-dwarfs or even brown dwarfs in the Galactic disk and bulge.

The formation theory of gas giants around the low-mass host remains to be fully elucidated. According to the core accretion theory, it is difficult to form gas-giant planets in the disks around low-mass stars (Pollack et al., 1996; Ida & Lin, 2004; Laughlin et al., 2004; Kennedy et al., 2006). Gravitational instability in the protoplanetary disk may play the important role in the formation of gas giants (Boss, 2006). In order to constrain the formation theory, more observational samples with low-mass hosts are required. By using the microlensing method, many planetary systems with low-mass stars have been discovered in orbital separation between $\sim 0.2-10$ au (Street et al., 2013; Skowron et al., 2015; Nagakane et al., 2017; Jung et al., 2018a). This is complementary to the other detection methods.

In this paper, we present the analysis of the planetary microlensing event MOA-bin-29 with a short Einstein radius crossing time of $t_{\rm E}\sim 4-7$ days, which suggests the host is a low-mass object. Section 2 explains observations and data. Our light-curve modeling method and result are shown in Section 3. In Section 4, we derive an angular Einstein radius from the source magnitude and color. In Section 5, physical parameters of the lens system are estimated with a Bayesian analysis. Finally, we discuss our analysis and reach conclusions in Section 6.

2 Observations

The Microlensing Observations in Astrophysics (MOA; Bond et al., 2001; Sumi et al., 2003) collaboration conducts a microlensing exoplanet survey toward the Galactic bulge by using the 1.8m MOA-II telescope with a 2.2 ${\rm deg^{2}}$ wide field-of-view (FOV) CCD-camera, MOA-cam3 (Sako et al., 2008) at Mt. John University Observatory in New Zealand. Thanks to the wide FOV, a high-cadence survey observation can be conducted. MOA survey uses a custom wide-band filter referred as $R_{\rm MOA}$, corresponding to a Cousins $R$- and $I$-band. The MOA photometry is reduced by using the MOA’s implementation of the Difference Image Analysis (DIA) pipeline (Bond et al., 2001).

The Optical Gravitational Lensing Experiment (OGLE; Udalski, 2003 ) also conducts a microlensing survey at the Las Campanas Observatory in Chile. The third phase of the survey, OGLE-III used the 1.3m Warsaw telescope with a 0.35 ${\rm deg^{2}}$ FOV CCD-camera. Currently, its forth phase, OGLE-IV (Udalski et al., 2015) started its high-cadence survey in 2010 with a 1.4 ${\rm deg^{2}}$ FOV CCD-camera. The OGLE photometry is reduced by the OGLE’s implementation of the DIA photometry pipeline (Udalski, 2003).

The gravitational microlensing event MOA-bin-29 reached to near peak on 2006 July 14 (${\rm HJD^{\prime}}$ $=$ ${\rm HJD}-2450000$ $\sim$ $3929.77$), at the J2000 equatorial coordinates $({\rm RA,Dec})$ = ($17^{h}$ $57^{m}$ $30^{s}$.$23,-29^{\circ}$ $44^{\prime}$ $11^{\prime\prime}.63$) corresponding to Galactic coordinates $(l,b)$ = $(0.633^{\circ},-2.636^{\circ})$. Figure 1 shows the MOA and OGLE-III light curves. This event was not detected in the real-time analysis but was found only after the off-line analysis of MOA database during $2006-2014$ (Sumi et al. 2019, in preparation). There are several possible reasons why this event was not detected by the MOA Alert system. First, this was a short duration event. Second, the alert system had just started since 2006 and the baseline was not long enough to distinguish from other variables. In this off-line analysis, the $2006-2014$ MOA Galactic bulge data have been re-analyzed and the events were detected using a criteria that is different from MOA alert system. Since this event was not alerted, there were no follow-up observations and only the survey data is available. Fortunately, however, the event is located in MOA field gb9, which was observed with the highest cadence of 10 minutes, so we have good coverage of this short event. The event is also located in the OGLE-III field BLG102, and we obtained data covering some part of the light curve during its magnification. Consequently, we must characterize the anomaly with only the survey groups’ observations. Figure 2 shows the reference image around MOA-bin-29. The green cross indicates the position of the event detected on the difference images. We found that the source star is much fainter than a nearby bright star (located up and to the left of the event location in Figure 2).

The photometric error bars produced by the data pipelines can be underestimated (or more rarely overestimated). We should consider other systematic errors caused by observational conditions and so on. In order to get proper errors of the parameters in the light-curve modeling, we empirically normalize the error bars by using the standard method of (Bennett et al., 2008). We use the formula,

[TABLE]

where $\sigma^{\prime}_{i}$ is the $i$ th renormalized error, $\sigma_{i}$ is the $i$ th error obtained from DIA, and $k$ and $e_{\rm min}$ are the renormalizing parameters. We set $e_{\rm min}=0.003$ to account for flat-fielding errors, and we adjust the value of $k$ and $\chi^{2}{\rm/dof}=1$. The normalization parameters and the number of data of each telescope are given in Table 1. We use the MOA data for $2006-2014$ and the OGLE data for $2006-2009$.

We investigated the possibility that this short magnification is not due to the microlensing, but the other artifacts or intrinsic variability of the star (Sumi et al., 2011; Bennett et al., 2012). First, we examined the pixel level DIA images of the target and confirmed that the event is not due to the fast moving objects nor cosmic ray hits. Second, we show the baseline of the full light curve except during the event at $3929.5<{\rm HJD^{\prime}}<3932.5$ in Figure 3. Here, we omitted the MOA data points with flux errors $>1000$ ADU or with seeing $>4$ pixels for clarity. We found that there are no other obvious magnifications in the baseline for $2002-2018$. This indicates that the event does not likely consist of a cataclysmic variables (CVs) nor flare stars because most of them repeat in a timescale of a few years or less. Additionally, the light-curve shape of these flare events is usually fast-rise and slow-decline, while the MOA-bin-29 light curve shows the opposite, i.e., slow-rise and fast-decline. Furthermore, the variables that can be modeled by binary microlensing light curves tend to have physically unlikely parameters and a small number of magnified data points (Bennett et al., 2012) and neither of these conditions apply to MOA-bin-29. These considerations indicate that this event is a microlensing event.

3 Light-curve Models

The single-lens light-curve model depends on three parameters: the time of lens-source closest approach $t_{0}$, the Einstein ring crossing time $t_{\rm E}$, and the impact parameter in units of the Einstein radius $u_{0}$. Binary-lens models require four additional parameters: the planet-host mass ratio, $q$, the planet-host separation in units of the Einstein radius, $s$, the angle between the trajectory of the source and the planet-host axis, $\alpha$, and the ratio of the angular source size to the angular Einstein radius. If we were to include microlensing parallax, we would need two additional parameters. The model flux $F(t)$ of magnified source as a function of time $t$ can be given by,

[TABLE]

where $A(t)$ is a magnification of the source flux at $t$, $F_{S}$ is the baseline flux of the source star, and $F_{b}$ is the baseline flux of any unresolved light. The flux on the DIA image, $\Delta$ $F(t)$, is the difference between $F(t)$ and the flux on the reference image at $t_{\rm ref}$, $F(t_{\rm ref})={\rm const.}$, and is given by

[TABLE]

The large number of the parameters of the microlensing event and correlations with each parameter make it difficult to search for the best-fit model parameters. The number of nonlinear parameters can be reduced by noting that equation 2 is linear in $F_{S}$ and $F_{B}$, so these parameters can be solved for directly for any $F(t)$ that is considered (Rhie et al., 1999).

The initial modeling of this event that led to its classification as a likely planetary microlensing event was done using three different modeling methods, using the modeling codes of Bennett (2010), Bozza (2010), and Sumi et al. (2010). Of particular importance for an event like MOA-bin-29, which does not have a significant signal from part of the light curve that resembles a single-lens light curve, is the variation of the Bennett (2010) method described in Bennett et al. (2012). This search code variation centers the grid search on the centers of caustics of wide or close binaries that are widely separated from the central caustic, and is able to efficiently find the correct solutions for events dominated by non-central caustics, such as MOA-bin-1 (dominated by a wide planetary caustic) and MOA-bin-3 (dominated by a minor image caustic of a close binary-lens system). This analysis found that only central caustic models were competitive within $\Delta\chi^{2}$ of 40.

After this preliminary analysis ruled out the widely separated caustic models, we used the Sumi et al. (2010) code for our detailed analysis. This code combines the Markov Chain Monte Carlo (MCMC) algorithm (Verde et al., 2003) with the image-centered ray-shooting method (Bennett & Rhie, 1996; Bennett, 2010). First, we performed a broad grid search over the at 9680 different grid points, over $(q,s,\alpha)$ space with other parameters free. Next, we refined all parameters for the best 100 models with the smallest $\chi^{2}$ to search for the global best-fit model. In conducting our grid search, we have set the initial parameters with range of $t_{0}\pm 10$ days, and $u_{0}\pm 0.5$. This, our grid search covers a wider parameter space. Furthermore, when we refine the model by freeing the $q$, $s$, and $\alpha$ parameters, the MCMC will probe for the wider solutions outside of each local minima. For example, the MCMC can find the solution with $q>1$ where $t_{0}$ becomes the time of the magnification near the secondary companion.

3.1 Limb Darkening

Binary-lens events usually have caustic crossings of cusp approaches that resolve the finite angular size of the source, so we must include the limb darkening of the source star. In order to take these effects into account, we adopt the following linear limb-darkening law:

[TABLE]

where $S_{\lambda}(\vartheta)$ is a limb-darkening surface brightness. The effective temperature of the source star estimated from the extinction-free source color presented in Section 4 is $T_{\rm eff}\sim 4939\,\rm K$ (González Hernández & Bonifacio, 2009). Assuming the surface gravity $\log g=4.5$ and metallicity of $\log[{\rm M/H}]=0$, we find limb-darkening coefficients of $u_{I}=0.5880$ and $u_{R}=0.6809$ from the ATLAS model (Claret & Bloemen, 2011). For the $R_{\rm MOA}$ passband, we use the coefficient for $u_{\rm Red}=0.6345$, which is the mean of $u_{I}$ and $u_{R}$.

3.2 Best-fit model and Degenerate models

By the grid search, we found the best binary-lens model and some local minima. The best binary model, wide1 model, is favored over the single-lens model by $\Delta\chi^{2}\sim 154$. Figure 1 shows a clear anomaly in the light curve from the single-lens model, and the best-fit binary-lens model can explain the anomaly near the peak. Figure 4 shows locations of these degenerate models in $q-s$ space from our five Markov chains of the $\chi^{2}$ distribution of the planet/host mass ratio $q$, and the planet-host separation, $s$, in the range of $4.0\times 10^{-3}\leq q\leq 1.0\times 10^{-1}$ and $0.4\leq s\leq 2.4$. The points are color coded based on $\Delta\chi^{2}$ from the best $\chi^{2}$ of the best model, the wide1 model. Thus, we took a closer look at these models. Figure 5 shows the caustic geometries of each model, and Figure 6 shows the light curves of all degenerate models. The parameters of models are listed in Table 2. Figures 7 shows zooms of the characteristic part of the light curves.

The five competing models are divided into wide models, with $s\ >\ 1$, and close models, with $s\ <\ 1$. The Einstein radii crossing times are short (t_{\rm{E}}$$\sim 4-7 days) for all models, and all models have a planetary mass ratio.

We now describe each of these degenerate models.

The wide1 model: This is the best-fit model with a planetary mass ratio of $q=1.6\times 10^{-2}$ and a separation of $s=1.2$. Figure 5 (a) and 6 shows the light curve and the caustic of this model. The Einstein radius crossing time $t_{\rm E}$ is only 4 days. We can see a bump around around $\rm HJD^{\prime}$ $\sim$ $3929.7-3930.0$ due to cusp crossing in Figure 7 (b) followed by the main peak due to the caustic exit where the MOA data have a good coverage in Figure 7 (a). Thanks to its caustic crossing feature, the clear finite source effects are detected by $\Delta\chi^{2}\sim 64$.

The wide2 model: Because the parameters such as $q$ and $s$ are slightly different from those of the wide1 model, the shape of the caustic is similar to that of the wide1 model in Figure 5 (b). We can see a bump similar to as that of the wide1 model on $\rm HJD^{\prime}$ $\sim$ $3929.7-3930.0$ because the source crosses the similar cusps as shown in Figure 7 (b). Although the source crosses a different part of the caustic from that of model 1, the features of the light curves of both models are alike during the data coverage.

The wide3 model: This model has a clearly different features in the caustic shape and the light curve from those of the wide1 model. This model is disfavored against the model 1 by $\Delta\chi^{2}\sim 7.1$. The mass ratio is $q=0.6\times 10^{-2}$ and the separation is $s=1.7$, which is larger than that of the wide1 model. Additionally, the Einstein radius crossing time is about 7 days, which is twice as long as that of the wide1 model. From the light curve in Figure 7 (c), we find another bump around $\rm HJD^{\prime}$ $\sim$ $3921.5-3923.5$ due to a cusp approach to a planetary caustic, and the main magnification arises by approaching to a cusp of the central caustic.

The close1 and close2 model: These two models have star–planet separations of $s<1$. The close1 and the close2 models are disfavored by $\Delta\chi^{2}\sim 13.2$ and $\Delta\chi^{2}\sim 13.9$, respectively. These two models have similar parameters, so the shapes of caustic geometry of both models are similar in Figure 5 (d) and (e). As for the close1 model, the mass ratio is $q=1.2\times 10^{-2}$ and the separation is $s=0.5$, which is smaller than that of the wide1 model. The Einstein radius crossing time is about 5 days. The light curves of both models are characterized by a cusp approach of the central caustic. The difference between the models is the source trajectory, and regarding close2 model, a bump around $\rm HJD^{\prime}$ $\sim$ $3936.8-3937.2$ comes from a planetary caustic crossing (Figure 7 (d)).

According to Figure 6 and 7, the difference between models are characterized by some bumps. The wide1 and the wide2 models have a bump at $\rm HJD^{\prime}$ $\sim$ $3929.7-3930.0$, and the wide3 model has a bump at $\rm HJD^{\prime}$ $\sim$ $3921.5-3923.5$. The close1 model does not have a bump, but the close2 model has a bump at $\rm HJD^{\prime}$ $\sim$ $3936.8-3937.2$.

The bump around $\rm HJD^{\prime}$ $\sim$ $3929.7-3930.0$ is most likely to be the real caustic feature for several reasons from the viewpoint of the data. First, the weather was clear at night during $\rm HJD^{\prime}$ $\sim$ $3929.80-3929.85$ and the seeing was also good. Around $\rm HJD^{\prime}$ $\sim$ $3929.80$, the flux ($\sim 3000$) was significantly larger than those of other small bumps ($\sim 1000$). Second, we checked the difference images during the bump and found that the difference images at $\rm HJD^{\prime}$ $\sim$ $3929.82$ indicate a variable object to be detected (Figure 8 (a)). The center of magnification of that image is consistent with those of the difference images around the main peak at $\rm HJD^{\prime}$ $\sim$ $3931$. Third, we found that the variable object in the difference image at $\rm HJD^{\prime}$ $\sim$ $3929.82$ is brighter than that of any in the different images during $\rm HJD^{\prime}$ $\sim$ $3929.82-3929.94$ (one of the images is in Figure 8 (b)). Therefore, the bump around $\rm HJD^{\prime}$ $\sim$ $3929.7-3930.0$ might have been caused by astrophysical origins.

However, it is difficult to judge if other small bumps (around $\rm HJD^{\prime}$ $\sim$ $3921.5-3923.5$ and $\rm HJD^{\prime}$ $\sim$ $3936.8-3937.2$) are due to the systematics on the baseline or the real caustic feature because these data points do not indicate any bad seeing, air masses, nor $\chi^{2}$ in PSF fitting. The $\chi^{2}$ differences due to these bumps are small, $\sim 10$ at most. We can see a lot of $2-3\sigma$ outliers that can cause similar bumps on the baseline data. The visual inspection of difference images during the bumps do not indicate any strange features such as saturation images, but a significant variable object cannot be detected in them. Although the effects of the airmass, the seeing, and the differential refraction for each data point have been corrected during the data reduction, it is not surprising that such small systematics remain in that dense stellar field.

Thus, we investigated whether the existence of these bumps affects the final results in the Appendix 6. As a result, we found that the final conclusion of the analysis remains the same whether we include these bumps or not, because they do not affect the MCMC distributions. However, the small $\chi^{2}$ differences cannot be used to compare the competing models. Instead, the range of parameters in these models, $6.0\times 10^{-3}\leq q\leq 2.0\times 10^{-2}$ and $0.5\leq s\leq 1.8$, should be taken as an uncertainty, conservatively. Therefore, we use these original five models without removing any the data points but recognize them as best-fit models equally, not weighting by $\Delta\chi^{2}$, in a Bayesian analysis in Section 5.

We also checked if the caustic feature producing the bump at $\rm HJD^{\prime}$ $\sim$ $3929.7-3930.0$ is real or not from the viewpoint of modeling, and reached two conclusions. First, we found new local minima with a similar bump around $\rm HJD^{\prime}$ $\sim$ $3929.7-3930.0$ when we conducted a new light-curve modeling grid search after removing all of the data points from that night. This implies that the bump is a real caustic feature. Second, this bump is closer to the center of mass, which is more likely to cross the caustics. On the other hand, the other bumps due to the planetary caustic are much farther from the center of mass, which means that the source is much more likely to miss these small caustics, particularly if planetary orbital motion or microlensing parallax was included in the modeling. So, the rest of the light-curve data implies that the bump at $\rm HJD^{\prime}$ $\sim$ $3929.7-3930.0$ is real, but it does not imply that the other light-curve bumps are real. Therefore, we conclude the bump at $\rm HJD^{\prime}$ $\sim$ $3929.7-3930.0$ is a real caustic feature.

We could detect clear signals of the finite source effects in the three wide models, but we obtained very weak measurements of $\rho$ for the close models. The best-fit close models are favored by only $\Delta\chi^{2}\sim 17$ and $\sim 6$, respectively, over models with $\rho=0$. We detect strong signals of the finite source effects for the other models.

3.3 Binary-source model

There is a possible degeneracy between the single-lens, binary-source model (1L2S) and the binary-lens, single-source model (2L1S) (Griest & Hu, 1993). For the 1L2S model, the magnification $A$ is expressed in the following equation:

[TABLE]

where $A_{1}$ and $A_{2}$ are the magnification of the two sources with flux $F_{1}$ and $F_{2}$, respectively, and $q_{F}$ is the flux ratio between the two sources ($=F_{2}/F_{1}$). For 1L2S model, the magnification $A$ depends on the wavelength unless the two source stars have the same color. By using the color difference expected for the two sources of unequal luminosity, the 2L1S/1L2S degeneracy could be solved (Gaudi, 1998). For some microlensing events, the 2L1S/1L2S degeneracy is broken with this method, and this confirms the planetary models (Shvartzvald et al., 2014; Zang et al., 2018). However, due to the poor data coverage for this event, the model parameters are uncertain, and this makes it difficult to use the color-shift method mentioned above to confirm the planetary interpretation more strongly. We searched for the best 1L2S model and found that the best 1L2S model is disfavored over 2L1S model by $\Delta\chi^{2}=\chi^{2}_{\rm 2L1S}-\chi^{2}_{\rm 1L2S}\sim 29$. Figure 9 shows the comparison between 2L1S/1L2S. Table 3 shows the parameters of the best-fit 1L2S model. According to Figure 9, the sizable fraction of $\chi^{2}$ differences arise from the bump, which we conclude is the real caustic feature. The $\chi^{2}$ difference from the outside of the bump is $\chi^{2}_{\rm 2L1S}-\chi^{2}_{\rm 1L2S}\sim 12.84$. So, the best-fit 1L2S model is disfavored over the 2L1S model both from the bump and from the rest of the light curve. Thus, we conclude that the possibility of the 1L2S model is excluded.

3.4 Microlensing parallax model

The microlensing parallax is an effect caused by the orbital motion of the Earth. Although it is known that there is little possibility of the detection of the parallax effect for such a short duration event (Gould & Loeb, 1992; Gaudi, 2012), we also considered a parallax model for completeness. Then, we found that the parallax model improves the fit only by $\Delta\chi^{2}\sim 8.40$, and the value of the parameters are $\pi_{\rm E,E}=407\pm 95$ and $\pi_{\rm E,N}=352\pm 88$, which are quite larger than the ordinary value ($<\ 1$). Therefore, we ruled out the parallax model.

4 Angular Einstein Radius

Thanks to the detection of the finite source effects, we can constrain the lens physical properties by estimating the angular Einstein radius $\theta_{\rm E}=\theta_{\rm*}/\rho$. We can get $\rho$ from the light-curve modeling and the angular source radius $\theta_{\rm*}$ by using empirical relation of $\theta_{\rm*}$, the intrinsic source color $(V-I)_{\rm S,0}$ and magnitude $I_{\rm S,0}$ (Boyajian et al., 2014). Because there is no $V$-band data during the magnification, we estimated the source color as follows. We cross-referenced stars in the MOA DoPHOT catalog of the reference image with the stars in OGLE-III photometry map (Szymański et al., 2011) within $120^{\prime\prime}$ around the source star. By using 91 cross-referenced stars, we derived the following color-color relation (see Figure 10):

[TABLE]

If we have a good measurement of the OGLE $I$-band source magnitude, $I_{\rm S,OGLE}$, we could derive the $(V-I)_{\rm S}$ from this formula. However, $I_{\rm S,OGLE}$ from the light-curve fitting has very large uncertainty because only a few data points during the low magnification are available.

Therefore, by following Bennett et al. (2008), we estimated the source color by taking the average color of main-sequence stars in Baade’s window observed by the $Hubble\ Space\ Telescope$, ($HST$, Holtzman et al. (1998)). In the color-magnitude diagram (CMD) (Figure 11), the black dots show the OGLE stars within $120^{\prime\prime}$ around the source star, and the green dots show the $HST$ stars that are adjusted for the reddening and extinction by using the the Red Clump Giants (RGC) color and magnitude, $(V-I,I)_{{\rm RGC,}HST}=(1.62,15.15)$ (Bennett et al., 2008). Because we do not have any good calibrated $I$-band source magnitudes, we derived its magnitude and color as follows. We solve for $I_{\rm S}$ and $(V-I)_{\rm S}$ using an iterative procedure. First, we estimate the initial source color, $(V-I)_{\rm S}$ from the average color of the main-sequence stars with the input $I_{\rm S}$ value. We then determine the new $I_{\rm S}$ values from this color and the $R_{\rm MOA}$ values from the light-curve model. After a few iterations, this converges. We used the $I_{\rm S,OGLE}$ value from the light-curve model for the initial $I_{\rm S}$ value.

We derive the extinction-free magnitude and color of source so as to calculate an angular Einstein radius by following a method similar to that of Yoo et al. (2004). The extinction-corrected magnitude can be determined from the magnitude and color of the centroid of the RGC feature in the CMD. In figure 11, the red point shows the centroid of RCG color and magnitude, $(V-I,I)_{\rm RCG}=(2.19,15.82)\pm(0.01,0.02)$ around the target. Assuming that the source star suffers the same reddening and extinction as the RGCs, we compare these values to the expected extinction-free RCG color and magnitude at this field of $(V-I,I)_{\rm RCG,0}=(1.06,14.41)\pm(0.07,0.04)$ (Bensby et al., 2013; Nataf et al., 2013), and as a result, we obtained the reddening and extinction by using the RCG color and magnitude, $(V-I,I)_{{\rm RGC,}HST}=(1.62,15.15)$ (Bennett et al., 2008).

We determined the angular Einstein radius for each model with the following method that we demonstrate with the parameters of the wide1 model. Assuming that the source star suffers the same reddening and extinction as the RGCs, we compare these values to the expected extinction-free RCG color and magnitude at this field of $(V-I,I)_{\rm RCG,0}=(1.06,14.41)\pm(0.07,0.04)$ (Bensby et al., 2013; Nataf et al., 2013), and as a result, we get the reddening and extinction to the source of $(E(V-I),A_{I})_{\rm RCG,0}=(1.13,1.40)\pm(0.07,0.05)$ for the best-fit model, wide1 model. Therefore, we determined the intrinsic source color and magnitude to be

[TABLE]

Then, we find the angular source radius with the empirical formula (Boyajian et al., 2014),

[TABLE]

where $\theta_{\rm LD}\equiv 2\theta_{*}$ is the limb-darkened stellar angular diameter, (Fukui et al., 2015). This relation is derived by using stars with colors corresponding to $3900\ <\ T_{\rm eff}\ <\ 7000$ (Bennett et al., 2017). We found the angular source radius $\theta_{\rm*}=0.45\pm 0.08\ {\rm\mu as}$ for the wide1 model, with the uncertainty dominated by the source color uncertainty rather than the $2\%$ uncertainly from the empirical formula.

Finally, we calculate the angular Einstein radius $\theta_{\rm E}=0.056\pm 0.012\ {\rm mas}$ and the lens-source relative proper motion $\mu_{\rm rel}=\theta_{\rm E}/t_{\rm E}=5.242\pm 1.144\ {\rm mas}\ {\rm yr}^{-1}$ for the wide1 model.

Since the source star color and magnitude depend on the model, the angular source radius also depends on the model, and we summarize the values of $I_{\rm S,0},(V-I,I)_{\rm S,0},\theta_{\rm*},\theta_{\rm E}$ and $\mu_{\rm rel}$ for each model in Table 4. As for the close models, we get only an upper limit for $\rho$, so we get only the lower limits of $\theta_{\rm E}$ and $\mu_{\rm rel}$.

5 Lens Physical Parameters by Bayesian Analysis

Because a microlensing parallax effect was not measured for this event, the lens mass cannot be directly measured from the light-curve models. In order to estimate the probability distribution of the lens properties, we conducted a Bayesian analysis (Beaulieu et al., 2006; Gould et al., 2006; Bennett et al., 2008) assuming the Galactic model of Han & Gould (1995) as a prior probability. We also assume a mass function used in Sumi et al. (2011) and extend it to the low-mass brown-dwarf regime ($0.001\leq M/M_{\odot}$). We use the measured $t_{\rm E}$ and $\theta_{\rm E}$ to constrain the lens physical parameters. The extinction-corrected blending flux, which includes the lens and unrelated ambient stars on the line of sight to the source star, is derived from the light-curve modeling and is set as the upper limit for the lens brightness. $H$-band magnitudes for the lenses of all models are estimated from the color-color relation of the main-sequence stars (Kenyon & Hartmann, 1995) and the isochrone model of 5 Gyr brown dwarfs (Baraffe et al., 2003). The extinctions are calculated from Cardelli et al. (1989).

Table 5 shows the summary of the lens physical parameters of each model, and we found that the results are divided into two types according to the mass of the host star. The median value of the probability distribution of the host mass for all models except wide3 indicates the host star is a brown dwarf, and that of the wide3 indicates that the host star is a late M-dwarf. As for the wide1 model, the lens system could be a gas giant with a mass of $M_{\rm p}=0.60\ M_{\rm Jup}$ orbiting a brown dwarf with a mass of $M_{\rm h}=0.03\ M_{\odot}$, located at $D_{\rm L}=7.12\ {\rm kpc}$ from the Earth, and a projected separation from the host star of $r_{\perp}=0.48\ {\rm au}$. If we assume a circular and randomly oriented orbit for the planet, the three-dimensional semi-major axis is expected to be $a_{\rm 3D}=0.59\ {\rm au}$. Assuming the wide3 model is correct, the lens system is likely a gas giant with a mass of $M_{\rm p}=0.63\ M_{\rm Jup}$ orbiting an M-dwarf with a mass of $M_{\rm h}=0.10\ M_{\odot}$, located at $D_{\rm L}=6.57$ kpc from the Earth, and a projected separation of $r_{\perp}=1.54\ {\rm au}$ and the three-dimensional semi-major axis is expected to be $a_{\rm 3D}=1.92\ {\rm au}$. Figure 12 shows the probability distribution of the lens parameters of the wide1 model and the wide3 model. The dark and light-blue regions show the $1\sigma$ and $2\sigma$ confidence intervals, respectively, and the vertical blue lines show the median value. The probability distributions of $V$-, $I$-, and $H$-band magnitudes with extinction of the host star of the wide1 model and the wide3 model are also shown in Figure 13.

Since the probability distributions of the lens physical parameters for all five models are consistent within $1\sigma$, we combined the distributions of these five models. Here, we combined these distributions with equal weight because all models are equally good within the uncertainty including the possible systematics. Figure 14 shows the combined probability distributions, which indicate that the lens system comprises a gas-giant planet with a mass of $M_{\rm p}=0.63^{+1.13}_{-0.39}\ M_{\rm Jup}$ orbiting a brown dwarf with mass of $M_{\rm h}=0.06^{+0.11}_{-0.04}\ M_{\odot}$ at $D_{\rm L}=6.89^{+1.19}_{-1.19}\ {\rm kpc}$ and a projected separation of $r_{\perp}=0.53^{+0.89}_{-0.18}\ {\rm au}$.

6 Discussion and Summary

We analyzed the short duration microlensing event MOA-bin-29, which was found only after conducting an off-line analysis of the MOA database using data from $2006-2014$ (Sumi et al. 2019, in preparation). Although we found five competing solutions, all degenerate models have a planetary mass ratio $\sim 10^{-2}$ and an Einstein radius crossing time of $~{}4-7$ days. The angular Einstein radius estimated from the detection of the finite source effects was used to constrain the lens parameters for some models. As a result of a Bayesian analysis, we found that the lens system is likely to be a gas giant orbiting a brown dwarf or a late M-dwarf.

Future high-resolution imaging with ground-based AO observations or space telescope could constrain the lens parameters (Bennett et al., 2006, 2007, 2015; Batista et al., 2014, 2015; Bhattacharya et al., 2017; Koshimoto et al., 2017; Bhattacharya et al., 2018).The source and the lens will be separated by $\sim 100$ mas for the wide1 model and $\sim 160$ mas for the wide3 model by 2025. According to the Figure 13 and Table 5, the lens $K$-band magnitude with extinction would be $K\sim 33$ mag and $K\sim 24$ mag for the wide1 and the wide3 model, respectively. These are $\sim 13$ mag and $\sim 4$ mag fainter compared to the source for the wide1 and the wide3, respectively. This indicates that there is only small chance to detect the lens flux even if the lens and sources are separated by $\sim 100$ mas. However the source flux could be detected using Keck AO (Batista et al., 2014) or $JWST$ (Gardner et al., 2006) because the measured source magnitude with extinction is bright enough to be detected. The source flux for each model is not well determined by the light-curve modeling, so the improvement of the accuracy of the source flux would constrain the degenerate models.

According to Suzuki et al. (2016), the detection efficiency and the survey sensitivity of planetary systems with a short duration is relatively low, so only a few microlensing planets with a short duration have been found (Bennett et al., 2014), and therefore any inferences drawn from these data have large statistical errors. It is therefore important to increase the number of samples of these planets. However the determination of lens properties for short duration events is difficult because in such a short event the measurement of a significant microlensing parallax effect is almost impossible. In order to solve the relation between lens physical parameters, other constraints are important, for example follow-up observations with high-resolution imaging as described above. Measuring the space parallax is also a powerful way to determine the lens parameters for short duration events (Shvartzvald et al., 2017; Chung et al., 2019). Especially, Earth-L2 separation between the ground and the $WFIRST$ telescope (Spergel et al., 2015) will be optimal for such short events with small Einstein radius.

The formation theory of a gas giant orbiting a brown dwarf and even that for a brown dwarf itself is still ambiguous. As a result of analyzing the event MOA-bin-29, we found the lens system is most likely to be a gas giant orbiting a brown dwarf, which is hard to form according to the core accretion theory. Therefore more and more accurate mass measurements of such systems are required. The microlensing method is a powerful way to detect these systems and some of such events have been discovered by the microlensing method (Han et al., 2013; Albrow et al., 2018; Jung et al., 2018a, b). Figure 15 shows the mass distribution of discovered systems with a brown dwarf/a late M-dwarf hosting a gas giant ($0.01\ M_{\odot}\ <\ M_{\rm h}\ <\ 0.3\ M_{\odot}$, $0.01\ M_{\rm Jup}\ <\ M_{\rm p}\ <\ 13.6\ M_{\rm Jup}$). For most of these events, the lens physical parameters such as mass are derived through a Bayesian analysis, which uses the mass function as a prior. The mass probability distribution of low-mass host stars ($M_{\rm h}\ <\ 0.1\ M_{\odot}$) depends strongly on the shape of the mass function which has large uncertainty. In Albrow et al. (2018), the lens parameters are directly derived by combining measurements from Earth and from the $Spitzer$ telescope (Gould et al., 2014). The increase in number of the mass measurements by the space parallax effect and high-resolution images are greatly anticipated by $WFIRST$ satellite in the near future, and would contribute significantly to clarification of the formation theory of a gas giant around a low-mass star.

T.S. acknowledges the financial support from the JSPS, JSPS23103002, JSPS24253004, and JSPS26247023. The MOA project is supported by the grant JSPS25103508 and 23340064. D.P.B., A.B., and D.S. were supported by NASA through grant NASA-NNX12AF54G. Work by N.K. is supported by JSPS KAKENHI grant No. JP15J01676. Work by A.F. is supported by JSPS KAKENHI grant No. JP17H02871. N.J.R. is a Royal Society of New Zealand Rutherford Discovery Fellow. A.S. is a University of Auckland Doctoral Scholar. The OGLE project has received funding from the National Science Centre, Poland, grant MAESTRO 2014/14/A/ST9/00121 to AU.

Appendix A Systematic Photometry Errors on the Baseline

In Section 3.3, our grid search found the five competing models. The $\chi^{2}$ differences between these different models are primarily due to the main light-curve peak, but several small bumps in the light-curve baseline also contribute. The bump at $\rm HJD^{\prime}$ $\sim$ $3929.7-3930.0$ is most likely to be a real caustic feature. The caustic geometries of the wide1 and wide2 models require a feature at approximately this time, and the significance of this feature is $\Delta\chi^{2}=15.1$. The other small bumps (Figure 7 (c) and (d)) could be due to systematic photometry errors. These are further from the peak, and the caustic geometry of the event does not require that the source trajectory will encounter them, particularly if realistic microlensing parallax and orbital motion are included. These features have a lower significance of $\Delta\chi^{2}=14.0$ and $\Delta\chi^{2}=7.3$ for features (c) and (d), respectively. On the other hand, these data are not associated with bad seeing, high airmass, or large PSF fit $\chi^{2}$ in the fit to the difference images. Visual inspection of the images does not show any obvious problems. We now investigate if the existence of these two bumps might affect the final results as follows. While we are not certain if these bumps are systematic errors or real features, we can estimate the possible uncertainty by this exercise.

In the best-fit models, we removed data points during night that contain the small bump causing the $\chi^{2}$ differences. Here, we removed all of the data points in the night because we cannot set a clear boarder line to cut these specific data points in question because there are a lot of similar $2-3\sigma$ outliers. Then, we conducted a grid search for the best-fit models. We repeated this process three times until there are no $\chi^{2}$ differences from the small bumps in the baseline. After this iteration, we found only two new local minima in addition to the models, which are almost same as the original best-fit five models. These new models are within the range of parameters of the original five models, $6.0\times 10^{-3}\leq q\leq 2.0\times 10^{-2}$ and $0.5\leq s\leq 1.8$. The distributions of the MCMC chains of these new local minima also roughly overlap to those of original models. Thus, the bumps that we decided are likely to be due to systematics influence only the best fit models. The best fit models are decided by very small changes in $\chi^{2}$. However, the inclusion, or not, of these bumps due to systematics is irrelevant because they do not affect the MCMC distributions, which show the real determination of the features of the light curve.

Furthermore, in order to check if the final probability distribution of the lens properties would change with these new local minima, we conducted a Bayesian analysis including these new models. Figure 16 shows the combined probability distributions which include the probability distribution of three new local minima and indicates the lens system comprise a gas-giant planet with a mass of $M_{\rm p}=0.61^{+1.12}_{-0.38}\ M_{\rm Jup}$ orbiting a brown dwarf with mass of $M_{\rm h}=0.05^{+0.11}_{-0.04}\ M_{\odot}$ at $D_{\rm L}=6.95^{+1.19}_{-1.16}\ {\rm kpc}$ and a projected separation of $r_{\perp}=0.51^{+1.05}_{-0.21}\ {\rm au}$. Thhese probability distributions of lens physical parameters are consistent with the original result in Section 5. We get a similar result when we conduct the same analysis, assuming the bump around $\rm HJD^{\prime}$ $\sim$ $3929.7-3930.0$ is due to the systematics, to be safe.

Therefore, we use these original five models without removing any data points but recognize the range of these models as uncertain. Thus, we treat these five models equally, not weighting by $\Delta\chi^{2}$, in a Bayesian analysis, in Section 5.