Blind extraction of an exoplanetary spectrum through Independent Component Analysis

1.2 MB·English

by  Ingo P. Waldmann

Blind extraction of an exoplanetary spectrum through
Independent Component Analysis
I. P. Waldmann, G. Tinetti
Department of Physics & Astronomy, University College London, Gower Street, WC1E 6BT, UK
3 [email protected]
1
0
2 P. Deroo

n Jet Propulsion Laboratory, California Institute of Technology, 4800 Oak Grove Drive, Pasadena, California
a 91109-8099, USA
J

7
1 M. D. J. Hollis, S. N. Yurchenko, J. Tennyson

Department of Physics & Astronomy, University College London, Gower Street, WC1E 6BT, UK
]
P
E
.
h
p
ABSTRACT
-
o
Blind-source separation techniques are used to extract the transmission spectrum of the hot-
r
t Jupiter HD189733b recorded by the Hubble/NICMOS instrument. Such a ‘blind’ analysis of
s
a the data is based on the concept of independent component analysis. The de-trending of Hub-
[ ble/NICMOS data using the sole assumption that nongaussian systematic noise is statistically

1 independent from the desired light-curve signals is presented. By not assuming any prior, nor
v auxiliary information but the data themselves, it is shown that spectroscopic errors only about
1 10 - 30% larger than parametric methods can be obtained for 11 spectral bins with bin sizes of
4
∼ 0.09µm. This represents a reasonable trade-off between a higher degree of objectivity for the
0
non-parametric methods and smaller standard errors for the parametric de-trending.
4
. Results are discussed in the light of previous analyses published in the literature. The fact
1
that three very different analysis techniques yield comparable spectra is a strong indication of
0
the stability of these results.
3
1
Subject headings: methods: data analysis — techniques: spectroscopic — planets and satellites: atmo-
:
v spheres — planets and satellites: individual(HD189733b)
i
X
1. Introduction 2010; Swain et al. 2008, 2009a,b; Tinetti et al.
r
a 2007b, 2010; Crouzet et al. 2012). The aim to
The field of exoplanetary spectroscopy is as
characterisesmallerandsmallerplanetsisequally
rapidly advancing as it is new. It has come
a quest for higher and higher precision measure-
from the first detection of spectroscopic features
ments, which are often limited by the systematic
in an exoplanetary atmosphere (Charbonneau et
noise associated with the instrument with which
al. 2002), to an ever more detailed characterisa-
thedataareobserved. Thisisparticularlytruefor
tionofavarietyoftargets. (e.g.Agoletal. 2010;
general, non-dedicated observatories. In the past,
Beaulieu et al. 2010, 2011; Berta et al. 2012;
parametric models have extensively been used by
Charbonneau et al. 2005, 2008; Deming et al.
several teams to remove instrument systematics
2005, 2007; Grillmair et al. 2008; Knutson et
(e.g. Agol et al. 2010; Beaulieu et al. 2010,
al. 2007a; Snellen et al. 2008, 2010a,b; Brogi
2011;Brownetal. 2001;Burkeetal. 2010;Char-
et al. 2012; Bean et al. 2011; Stevenson et al.
1
bonneau et al. 2005, 2008; Deming et al. 2007; ment and stellar noise components. ICA belongs
D´esert et al. 2011; Gibson et al. 2011; Grillmair to the category of unsupervised machine learn-
et al. 2008; Knutson et al. 2007a; Pont et al. ing algorithms (i.e. the algorithm does not need
2008;Singetal. 2011;Swainetal. 2008,2009b). to be trained prior to use). It does not require
Parametric models approximate systematic noise auxiliary information such as OSVs but only the
via the use of auxiliary information of the instru- observed data themselves. This approach is also
ment, the so called optical state vectors (OSVs). known as blind decorrelation as the assumptions
Such OSVs often include the X and Y-positional on the system are minimal. As described later
drifts of the star or the spectrum on the detector, on, ICA assumes a linear combination of inde-
the focus and the detector temperature changes, pendent components (ICs) to form its systematic
aswellaspositionalanglesofthetelescopeonthe noise model.
sky. ByfittingalinearcombinationofOSVstothe The issue of poorly constrained parameter
data, the parametric approach derives its system- spaces is in fact not new in astrophysics and has
atic noise model. We refer to this as the ‘linear, given rise to an increased interest in blind-source
parametric’method. Inthecaseofdedicatedmis- separation algorithms. Cosmological and extra-
sions,suchasKepler(Boruckietal. 1996;Jenkins galacticobservations,inparticular,areoftenanal-
et al. 2010), the instrument response functions ysedthroughfullyblind, non-parametricmethods
are well characterised in advance and conceived and ICA has successfully been used to separate
to reach the required 10−4 to 10−5 photometric the cosmic microwave background (CMB) or the
precision. For general purpose instruments, not signatures of distant galaxies from their galactic
calibrated to reach this required precision, poorly foregrounds (e.g. Chapman et al. 2012; Stivoli
sampled optical state vectors or a missing param- et al. 2006; Maino et al. 2002, 2007; Wang et
eterisationoftheinstrumentoftenbecomecritical al. 2010). Aumont & Mac´ıas-P´erez (2007), for
issues. Eveniftheparameterisationissufficient,it instance,canseparatetheinstrumentalnoisefrom
is often difficult to determine which combination the desired astrophysical signal by using ICA.
of these OSVs may best capture the systematic
In W12 the ICA approach was tested by pre-
effects of the instrument.
senting two single, de-correlated light-curves of
Given the intricacies of a parametric approach, the primary transit of HD189733b and XO1b ob-
several groups have worked towards alternative tained with Hubble/NICMOS in its spectroscopy
methods to decorrelate the data from instrumen- setting. In this paper, we apply the ICA de-
tal and stellar noise. The issue of time-correlated trending method to the extraction of an exoplan-
systematics in exoplanetary time series was dis- etary spectrum; we then compare the results ob-
cussed by Pont et al. (2006). Carter & Winn tainedtothepreviousonespublishedinthelitera-
(2009)developedawaveletbased,non-parametric, ture and discuss the advantages and limits of this
de-trendingof1/f noisecontaminatedlightcurves. technique for analysing exoplanetary spectra and
Thatte et al. (2010) proposed a selective prin- light-curves.
cipal component filtering to reduce instrument
and telluric systematics. More recently, Gibson 1.1. TheHD189733bHubble/NICMOSdata-
et al. (2012),here after G12, presented a non- set
parametric de-trending approach based on Gaus-
As explained in the introduction, the main fo-
sianProcesses(GP;Rasmussen&Williams 2006).
cus of our paper is the application and critical
GP, as implemented by G12, belongs to the class
discussion of the ICA method to de-trend time-
ofnon-parametric,supervisedmachinelearningal-
resolved spectroscopic data. We chose the transit
gorithms. It uses the observed data and OSVs
spectrum of the hot-Jupiter HD189733b recorded
to derive an optimal, non-linear systematic noise
byHubble/NICMOSasatestbed,themainreason
model for the observed data.
beingthatthesedatahavebeenanalysedbyother
In Waldmann (2012), here after W12, we
teams in the past and the results of the analysis
proposed independent component analysis (ICA;
have been quite debated in the literature.
Hyv¨arinen 1999) as an effective way to de-
The NICMOS data were first published by S08
correlate the exoplanetary signal from the instru-
2
and then re-analysed by Gibson et al. (2011) ual source signals as time series, i.e. k = l.
using a similar linear parametric de-trending ap- Assuming only one source signal is astrophysi-
proach. Gibson et al. (2011), failing to retrieve cal and one is Gaussian per time series, we can
the transmission spectrum reported by S08, at- state the total number of source signals to be
tributed the discrepancy to a high degree of de- N = N + 2, where N is the number of
sn sn
generacy between the systematic noise correction systematic noise source signals.
model used and the extracted spectrum. Swain Itisworthnotingthatforovercompletesets,we
et al. (2011) argued, though, this discrepancy havemoretimeseriesavailablethansourcesignals
being due to poorly derived OSVs by Gibson et contained within, i.e. k > l. In these cases we
al. (2011). A more recent re-analysis of the data can reduce the dimensionality of the data set (for
byG12usingnon-parametric,non-linearGaussian example using principal component analysis) or
Processes,foundasolutionconsistentwiththatof select a sub-set of estimated source signals given
S08, but error bars on average 2.4 times as large. some selection criteria (W12).
Given these debates, it becomes critical to un- We can also express eq. 1 in matrix form as:
derstandhowfarwecanpushtheanalysiswhenno
prior or auxiliary information for the instrument
x=As (3)
is assumed.
where x is the column vector containing the mea-
2. Data analysis sured time series, x , i.e. x = (x ,x ,...,x )T,
k 1 2 N
s is the column vector of independent source sig-
In this section we briefly introduce the ICA
nals, s = (s ,s ,s ,s ,...,s )T. We may
a wn sn1 sn2 N
technique and we apply it to the primary tran-
also write s = s +s +s to clearly differen-
a wn sn
sit data of HD189733b, as recorded by Hub-
tiate between the astrophysical, white noise and
ble/NICMOS.
systematic components. A is the N ×N dimen-
sional ‘mixing matrix’ comprised of the weights
2.1. Independent Component Analysis
a .
k,l
Let us assume we have multiple, simultaneous ThemotivationofICAistoestimateAwithout
observations of a mixture of signals. In time- priorknowledgeofAors(Hyv¨arinen 1999). This
resolvedspectroscopy(spectrophotometry)weob- is achieved by making the stringent assumption
tainalight-curvesignaloftheexoplanetarytransit thatthesourcesignalscomposings,s ,arestatis-
k
per resolution element of the spectrograph. Here, tically independent from each other. Such an as-
each measured time series is denoted by xk(t), sumptionisvalid,asoneexpectstheastrophysical
where k is its index in the individual time series light-curvesignaltobeindependentinoriginfrom
and N is the total number of time series. The systematic instrumental noise. ICA algorithms
measured signal can be assumed to be the sum hence attempt to de-compose observed signals, x,
of the astrophysical light-curve signal, sa(t), the to a set of independent source components, s. To
instrumental and stellar noise sources, ssn(t) and maximise said independence, several approaches
the white noise, swn(t). So: have been proposed (for a comprehensive sum-
mary: Hyv¨arinen & Oja. 2001; Comon & Jutten
2000). HerewefollowW12anduseavariantofthe
x (t)=a s (t)+s (t)+a s (t)+ (1)
k k,1 a wn k,2 sn1 FastICA algorithm, which maximises statistical
+ak,3ssn2(t)+...+ak,Nsnssn(t) independenceoftheestimatedsourcesignalssl by
maximising the nongaussianities of their respec-
or as sum of vectors (the time-dependance has
tive probability distributions (Hyv¨arinen 1999;
been dropped for clarity):
Hyv¨arinen et al. 1999; Hyv¨arinen & Oja 2000;
N(cid:88)sn Koldovsky´ et al. 2006; Hyv¨arinen & Oja. 2001;
xk =ak,1sa+ak,2swn+ ak,lssn (2) Comon & Jutten 2000).
l=3 A linear combination of individual source sig-
where l is the estimated source signal index. For nals, as in eqs. 1 and 2, is hence a valid assump-
non-overcomplete sets, we have as many individ- tion as each component is assumed to be fully
3
independent. However, a few subtleties regard- low pass filters, wavelet based approaches, etc.)
ing this approach are worth mentioning. Only butwedecidedtointerfereaslittleaspossiblewith
the signals that are common to all time series, the original data and opted for a slightly coarser
x (t), can realistically be de-convolved. Taking binning instead. This resulted in 11 light curves
k
Hubble/NICMOS as example, systematics intro- across the G206 grism band. We found the first
duced by the grism are easily detectable, as these of the 5 orbits to be very noisy and negatively
are ‘common’ to all dispersed wavelengths. Sim- impacting the efficiency of the algorithm and ex-
ilarly, detector-wide flat-fielding gradients are de- cludedthefirstorbitfromallfurtheranalysis. An
trendable. Ontheotherhand,localisedinter-pixel exampleofthe‘raw’light-curves’quality,at∼2.33
variationsarenotrepresentedinalltimeseriesand µm, can be found in figure 1.
wouldnotbede-correlatedwithICA.Theseprop- As described in W12, we used the extracted
erties lead to the discussion of ‘global’ and ‘local’ light-curves as input to the ICA algorithm to cal-
noise models later on in the text. Whilst the ef- culatethemixingmatrix,A,andits(pseudo)inverse,
fect of limb-darkening is reduced in the infra-red, the de-mixing matrix W = A−1. Once the de-
it is true that wavelength varying limb-darkening mixing matrix had been determined, the algo-
coefficients can impair the directextraction of the rithm tested the estimated components for their
lightcurve feature. However, as described in the nongaussianity and returned four main system-
nextsection,thesystematicnoisemodelmakesno aticnoisecomponentswhichdonotcorrelatewith
directuseoftheastrophysicalcomponent,sa,and the expected light-curve morphology. These com-
hence circumvents this potential limitation. ponents, comprising s , were extracted over the
sn
entire spectral range of the grism (referred to as
2.2. Application to HD189733b
‘global’ below) and showed a good degree of sep-
aration (figure 2, left side). ICA estimates the
The primary transit of HD189733b was ob-
mixing matrix A up to a sign and scaling factor,
served using HST/NICMOS in the G206 grism
meaning the source signals, s, are recovered but
setting and spanning five consecutive orbits. The
lack an overall scaling constant. In an analogy to
selected grism covers the spectral range of ∼1.51
principalcomponentanalysis,wecanthinkofthis
- 2.43 µm, see S08. The HST-pipeline calibrated
data were downloaded from the MAST1 archive as recovering the eigenvectors but not the eigen-
values. Due to this ambiguity, we need to deter-
and the spectrum was extracted using both stan-
dardIRAF2 routinesaswellasacustombuiltrou- minethescalingoftheindividualsystematicnoise
components,pertimeseriesx ,separately. Thisis
tine for optimal spectral extraction. Both extrac- k
done by fitting a systematic-noise-model (SNM),
tions are in good accord with each other but the
m , to the out-of-transit data of each time series
custom built routine was found to yield a better k
x , where m is the sum of all scaled systematic
signal to noise and was subsequently adopted for k k
noisecomponentsins ,i.e. m =(cid:80)Nsno s ,
all further analysis. In order to minimise inter- sn k l=1 k,l sn,l
where o is the scaling factor for the systematic
and intra-pixel variability of the NICMOS detec- k,l
noise signal s for a given time series k. A
tor, the instrument was slightly de-focused to a sn,l
Nelder-Mead minimisation algorithm (Press et al.
full-width-half-maximum(FWHM)of∼5spectral
2007) was used to fit for o . The scaling am-
channels per resolution element. This sets a limit k,l
plitude of each component for each light-curve is
on the maximum resolution, R, achievable.
given in figure 2 (right side). Once m is deter-
ICAislimitedbytheamountofGaussiannoise k
mined,wesubtractitfromtherawdatatogetthe
present in the data (Hyv¨arinen & Oja 2000). We
corrected time series y = x −m , see figure 1
found the minimal binning of 5 channels to be k k k
for an example.
too low in SNR and used a binning of 8 spectral
channels (∼ 0.09µm). Several data pre-processing
2.3. Lightcurve fitting and Error-bars
methods exist to decrease the Gaussian compo-
nent of time series data (e.g. kernel smoothing, Havingobtainedthede-trendedtimeseries,y ,
k
we proceed to model-fit these using the analytical
1http://archive.stsci.edu/ modelbyMandel&Agol (2002),fromhereMA02,
2http://iraf.noao.edu/ withorbitalandlimb-darkeningparameterstaken
4
from S08. The S08 quadratic limb-darkening pa-
rameterareinterpolatedtothecoarserwavelength
grid of this analysis. The transit depth, δ , is left
k
asonlyfreeparameter. Thetransitdepthsofall11
light-curvesconstitutetheexoplanetaryspectrum.
The model-fitting was performed using a Markov
Chain Monte Carlo (MCMC) algorithm and cross
checked using two variants of a Bootstrap Monte
Carlo analysis (see appendix A).
2.3.1. Markov Chain Monte Carlo
MCMC (Press et al. 2007) has become the
standard fitting routine for exoplanetary time se-
ries and radial velocity data (e.g. Ford 2006;
Burke et al. 2007; Bakos et al. 2007; Knutson
et al. 2007a; Cameron et al. 2007; Charbonneau
Fig. 1.— Raw light-curve at ∼2.33 µm (black et al. 2009; Bean et al. 2011; Kipping & Bakos
crosses),itsrespectivesystematicnoisemodel(red 2011;Gregory 2011;Crouzetetal. 2012;Knutson
squares), mk(t), composed out of the systematic et al. 2012). In this analysis we use an adaptive
componentsinfigure2. Thede-trendedfinallight- version of the standard Metropolis-Hastings algo-
curve is shown underneath (blue circles) with a rithm (Haario et al. 2001, 2006; Metropolis et
Mandel & Agol (2002) fit overlaid. al. 1953; Hastings 1970; Press et al. 2007). The
onlyfreeparameter, δ , wassettohaveauniform
k
prior ranging from R /R = 0 - 1. As we only
p ∗
fit for the transit depth we are not concerned by
inter-parameter correlations in this analysis.
WeperformaninitialMA02modelfit,m (t),
MA02,k
using a Nelder-Mead minimisation algorithm
(Press et al. 2007) and calculate the model sub-
tracted residual, r (t) = x (t)−m (t). We
k k MA02,k
take the fitted transit depth as starting value of
theMCMCchainandcalculatethevarianceofthe
normal sampling distribution as
(cid:88)
Fig. 2.— LEFT: Four retrieved nongaussian sys- σk2,sample =var(rk(t))+2 cov[rk(t),rk(t+τ)]
tematic noise components in the order of impor- τ=1
(4)
tance. They were computed over the whole spec-
tralrangeoftheG206grismanddescribethesys- where var is the variance of the residual and cov
tematic noise (instrumental and/or stellar) com- the auto-covariance for a given lag τ. This ac-
montoallspectralchannels. RIGHT:Scalingfac- counts for remaining autocorrelated noise in the
tors, o , of the systematic noise components on time series data. The MCMC algorithm was con-
k,l
the left. The colour coding is identical for both sequently run for 2×104 iterations to guarantee
plots. We can see that the first component at a good coverage of the posterior distribution, of
2.06 µm is sharply deviating from its own pass- whichexamplesareshowninfigures6a&7a. Here
band mean and the mean of all components at the error bars are the standard deviation of the
2.06 µm. This can indicate that the ‘global’ sys- posterior distribution.
tematicnoisemodeldoesnotdescribewellthesys- In addition to the MCMC algorithm described
tematics in this channel and could over-correct. here, we also estimated the standard error us-
ing two variants of a Bootstrap Monte-Carlo al-
5
gorithm, see appendix A. We find the retrieved
transit-depthsanderrorstobeingoodagreement
E[(cid:80)N G2 ]
with the MCMC method. σ = l=1,l(cid:54)=k kl , k,l=1,2,...,N. (7)
l E[G2 ]
kk
2.3.2. Source signal separation error
However, the original mixing matrix, A, and
The two core algorithms used in this analysis
the original source signals, s , are not generally
EFICA (Koldovsky´ et al. 2006) and WASOBI k
known for real data sets and equation 7 is only
(Yeredor 2000; Tichavsky´ et al. 2006b) can be
useful in the case of simulations. Tichavsky´ et al.
shown to be asymptotically efficient, i.e. reaches
(2006a); Koldovsky´ et al. (2006) and Tichavsky´
theCramer-RaoLowerBound(CRLB)inanideal
et al. (2008) have shown that despite the origi-
casewherethenonlinearityG(.)equalsthesignal’s
nal mixing matrix being unknown, an asymptotic
scorefunction. Inotherwords,thealgorithmsem-
estimate of the SIR can be made. A derivation
ployed here can be shown to converge to the cor-
of this process is beyond the scope of this paper
rect solution given the original source signals and
and we refer the interested reader to the relevant
in the limit of N →∞ iterations.
iter literature.
In reality the number of iterations is finite and
We now have the retrieval error on the individ-
imperfect convergence results in traces of other
ual source components, σ . In order to obtain the
sources to remain in the individual signals com- l
overall error on the systematic noise model, m ,
prising s. We can hence state that the estimated k
calculatedinsection2.2,wecomputetheweighted
de-mixing matrix, W is only approximately equal
sumoftheindividualsourceseparationerrorswith
totheinverseoftheoriginalmixingmatrix,A,i.e.
thepreviouslyretrievedsourcecomponentweight-
ing factors, o
W(cid:39)A−1 (5) k,l
This requires us to calculate the signal separation (cid:32)N(cid:88)sn (cid:33)1/Nsn
σ = o2 σ2 (8)
error(SSE)oftheanalysis. Ameasureofthiserror k,SSE k,l l
is the deviation of WA from the unity matrix by l=1
inspectingthevarianceofitselements(Koldovsky´
2.3.3. SNM fitting error
et al. 2006; Hyv¨arinen & Oja. 2001).
In addition to the above determined errors, we
To assert a good degree of separation, we can
also include an ICA fitting error which accounts
define G as the gain matrix. For a perfectly es-
for possible over-corrections of the global SNM
timated de-mixing matrix, W, the gain matrix is
to individual, poorer constraint light-curves. This
equal to its identity matrix
term becomes non-zero when the scaling of a sys-
tematicnoisecomponent,o (figure 2,righthand
G=WA=I (6) k,l
side), shows a 3σ significant deviation from the
In signal processing, the performance of blind- meanscalingofallotherlight-curves,o¯k. Inother
source separation algorithms is usually measured words, we expect the scaling of an individual sys-
bythesignaloverinferenceratio,SIR3.TheSIRis tematicnoisecomponent,ssn,tobeaslowlyvary-
the standard measure in signal processing of how ing function over wavelength for ‘globally’ esti-
well a given signal has been transmitted or de- mated systematic-noise components. If individual
convolved from a mixture of signals. Given the light-curves show a significantly larger positive or
inference to be a noise source, we can equate this negative scaling than expected, for an individual
to the more commonly used signal-to-noise ratio lightcurve, we can assume the nongausian noise
(SNR).Wecannowcalculatetheseparationerror of the affected light-curve not to be properly ac-
of the estimated source signal, s , in relation to counted for by this global model. In larger data-
l
the original source signal, s , using sets it is easier to exclude the affected light-curve
k
fromanyfurtheranalysis,whilstinsmalldata-sets
3standardliteraturealsoreferstoitsinverse,theinterference wetaketheamplitudeofthescalingfromitsmean
oversignalratio,ISR scaling as the error, i.e.
6
σ =|(o −o¯ )| (9)
k,SNM k,l k
In this analysis we find the ICA fitting error to
bezeroforallwavelengthsbutthe2.06µmspectral
point.
2.3.4. Final Error Bar
In summary the final error bar per time series,
k, consists of:
1. Standard error: Estimating the variance
inretrievedtransitdepthwhenmodelfitting
Fig. 5.— Comparision between ICA derived non-
of the de-trended lightcurves.
gaussian components on the left with convention-
2. Signal separation error: Estimatingthe allyderivedOSVsasfoundbySwainetal. (2008).
ICA component separation error.
3. Systematic noise model error: Esti- noise component is indicative of an overcorrection
mating errors due to noise model over- or which is reflected in the error-bar as described in
under-fitting for individual time series. the previous section. Overall, the error-bars re-
ported here are ∼ 10 - 30% larger than those re-
We now define the final error to be the sum of ported by S08 but ∼ 10 - 50% smaller than those
squares of the above mentioned error sources: reported by G12. It should be noted that this
analysisusesaslightlycoarserbinsizeyielding11
data points whilst the S08 and G12 analyses yield
(cid:113)
σ = σ2 +σ2 +σ2 18 spectral points.
k,TOTAL k,MCMC k,SSE k,SNM
(10)
4. Discussion
3. Results In this analysis we have computed the global
SNMandfoundagoodagreementwithpreviously
Figure1showstherawlight-curveat∼2.33µm
published results. Figure 4 shows the compari-
in black crosses with its corrected counterpart
son between the spectrum derived in this analy-
(blue circles) offset below. Here, much of the au-
sis, compared to the linear, fully parametric ap-
toregressive noise in the original data could be
proachbyS08andthenon-linear,non-parametric,
captured by the noise model (red squares) and re-
Gaussian Processes, approach by G12. Both non-
movedfromthefinalresult. Allrawandcorrected
parametric approaches yield slightly higher error
light-curves are shown in figure 3 with their re-
bars than the parametric approach. We find the
spective systematic noise models. The resulting
error on the ICA derived spectrum to be ∼10 -
spectrum is presented in figure 4 and table 1. We
30% bigger than those reported by S08 using a
find the retrieved spectrum to be in good agree-
coarser bin size of ∼ 0.09µm. These differences
ment with the ‘parametric’ analysis S08 and G12
in error bars between the ‘blind’ and ‘informed’
in terms of spectral shape which show-cases the
(meaning the use of auxiliary information of the
robustness of this methodology and the stability
instrument) approaches are not surprising. By
of the result as a whole.
not assuming any knowledge of the data or in-
The underlying noise of the spectral point at
strument, we are actively neglecting auxiliary in-
∼2.06 µm was flagged by the algorithm to be dis-
formation helpful to the de-trending of the data
crepant with the global systematic noise model.
set.
This can also be seen by the higher systematic
With the linear, non-parametric blind ICA
noise remaining in the corrected lightcurve (5th
method presented here, the uncertainties grow by
fromthetopinfigure3). Herethefirstsystematic
7
Fig. 3.— Left: Raw light-curves from 1.51µm (bottom) to 2.43µm (top) with fitted Mandel & Agol (2002)
model overlaid. The lightcurves were off-set for clarity. Middle: Systematic noise model for each raw
lightcurve in the left panel. Right: Final de-trended light-curves with fitted Mandel & Agol (2002) model
overlaid.
Fig. 4.— Final spectrum (red circles) obtained with the ICA algorithm described here and in Waldmann
(2012), overlaid on the ersults of Swain et al. (2008) (squares) and G12 (triangles).
8
Table 1
NICMOS transmission spectrum of HD189733b for a ‘global’ systematic noise model
correction and plotted in figure 4. The columns are wavelength (λ), planet-star-ratio
(R /R ) and the respective error-bar (σ(R /R )).
p s p s
λ (µm) R /R σ(R /R )×10−4
p s p s
2.429 0.155187 3.69383
2.336 0.155125 3.03754
2.244 0.155437 2.83050
2.152 0.155000 3.01554
2.060 0.155200 12.29019
1.968 0.155187 2.94831
1.876 0.155375 3.09120
1.784 0.155187 3.05973
1.691 0.154250 2.54556
1.599 0.154000 3.75548
1.507 0.156500 3.28963
upto30∼40%comparedtothelinear,parametric In other words, at least as many observed time
analysisbyS08(accountingforthelargerbinsizes series, x , are required as input to the algorithm
k
inthisanalysis)andafurther∼70%whenfurther than total number of nongaussian components in
relaxing the linear assumptions as for the non- the data. In this case, the minimum ‘local’ SNM
linear, non-parametric Gaussian Processes (G12). would include 5 spectral bins (4 systematic noise
In other words, we are trading smaller error-bars and1astrophysicalcomponent). With11spectral
for a higher degree of objectivity. bins in total, the NICMOS data-set is too small
In figure 5 we show a comparison of the ICA for this approach. However, for larger sets, this
derived nongaussian systematic noise components two stage ‘global’ + ‘local’ detrending becomes
and parametric OSVs. S08 identified the X and a viable solution. Furthermore, multiple obser-
Y-positional drifts on the detector to constitute vations of a transit/eclipse event with the same
the most important OSVs in their de-correlation instrumental setup can be helpful to further de-
whilst other OSVs are less significant. Similar correlate the observed data. We expect in fact
independent components are also present in this the astrophysical signal being stable throughout
analysis,whilsttheotherindependentcomponents consecutive transits, while the instrumental noise
differ more significantly from the parametric ap- being largely uncorrelated between observations
proach. The agreement between the results of (Waldmann et al. 2012).
this analysis and those available in the literature,
showcase the stability of this exoplanetary spec- 5. Conclusion
trum.
Here we present a reanalysis of a HD189733b
We find the global SNM approach to be most
primary eclipse spectrum from ∼1.51 - 2.43 µm
sensitive to slowly varying systematic trends
obtainedwithHubble/NICMOS.Thisanalysisdif-
across the data-set while local nongaussian de-
fersfrompreviouspublicationsinthatitusesblind
viations tend not to be captured. It is therefore
machine-learning to de-trend data with no prior
possible to generate, additionally, a ‘local’ SNM
or auxiliary information about the data or in-
for a sub-set of spectral bins or before binning
strument. Such blind-source de-convolution algo-
on the individual ‘raw’ spectral channels, should
rithms can be used alone or in conjunction with
the SNR permit it. It is important to remember
otherde-correlationtechniques,makingthemvery
that k ≥ N +1 as the input to the algorithm.
sn powerfultoolsintheanalysisofexoplanetarydata.
9
Thisisespeciallytrueforinstrumentsthatlackan Baluev, R. V., 2009, MNRAS, 393, 969
a priori calibration plan at the level of 10−4 pho-
Bakos, G. ., Noyes, R. W., Kovcs, G., et al., 2007,
tometric precision, as needed for this field.
ApJ, 656, 552
We compare our results with previously pub-
lishedanalysesofthesamedataset: anothernon- Barber, R. J., Tennyson, J., Harris, G. J., &
parametric, non-linear approach (G12) and the Tolchenov, R. N. 2006, MNRAS, 368, 1087
classical linear parametric method (Swain et al.
2008). We find that the error-bars of this analy- Bean, J. L., Miller-Ricci Kempton, E., Homeier,
sis are ∼ 10 - 30% larger than those reported by D., 2011, Nature, 478, 41
Swain et al. (2008). We attribute this difference
Bean,J.L.,D´esert,J.-M.,Kabath,P.,etal.,2011,
to the higher amount of auxiliary information in-
ApJ, 743, 92
jected in the parametric approach. Ultimately, it
is a trade-off between a higher degree of objectiv- Beaulieu, J. P., Kipping, D. M., Batista, et al.,
ity for the non-parametric methods and smaller 2010, MNRAS, 409, 963
errors for the parametric de-trending. Additional
Beaulieu, J.-P., Tinetti, G., Kipping, D. M.,
observations would have allowed much smaller er-
Ribas, et al., 2011, ApJ, 731, 16
ror bars and a more robust determination of the
signal at ∼2.06 µm.
Berta, Z. K., Charbonneau, D., Dsert, J.-M., et
Thefactthatthreeverydifferentanalysistech- al., ApJ, 747, 35
niquesyieldcomparablespectraisastrongindica-
tionofthestabilityoftheseresults. Theerrorbars Borysow, A., 2002, A&A, 390, 779
estimated in this paper through ICA are smaller
Borucki, W. J., Dunhm, E. W., Koch, D. G.,
than the ones of Gibson et al. (2012) through
Cochran, W. D., et al., 1996, Astrophys. &
GaussianProcesses,suggestingmoreinvestigation
Space Science 241, 111
is needed to identify the most effective techniques
to de-trend exoplanet atmosphere data and, most Brogi M., Snellen, I. A. .G., de Kok, R. J., et al.,
importantly, understand their limitations. 2012, Nature, 486, 502
Brown, T. M., Charbonneau, D., Gilliland, R. L.,
The authors thank Mark Swain and Filipe Ab-
Noyes, R. W. & Burrows, A., 2001, ApJ, 552,
dalla for helpful discussions. This work is sup-
699
ported by ERC Advanced Investigator Project
267219,STFC,NERC,UKSA,UCLandtheRoyal Burke, C. J., McCullough, P. R., Valenti, J. A.,
Society. All of the data presented in this paper Johns-Krull, C. M., Janes, K. A. et al., 2007,
were obtained from the Multimission Archive at ApJ, 671, 2115
the Space Telescope Science Institute (MAST).
STScI is operated by the Association of Universi- Burke, C. J., McCullough, P. R., Bergeron, L. E.,
tiesforResearchinAstronomy,Inc.,underNASA et al., 2010, ApJ, 719,1796
contractNAS5-26555. SupportforMASTfornon-
Chapman, E., Abdalla, F. B., Harker, G., et al.,
HSTdataisprovidedbytheNASAOfficeofSpace
2012, MNRAS, 423, 2518
Science via grant NNX09AF08G and by other
grants and contracts. Carter J.A., Winn J.N., 2009, ApJ, 704, 51
REFERENCES Collier Cameron, A., Wilson, D. M., West, R. G.,
Hebb, L., et al., 2007, MNRAS, 380, 1230
Agol, E., Cowan, N. B., Knutson, H. A., Deming,
D., et al., 2010, ApJ, 721, 1861 Charbonneau, D., Brown, T. M., Noyes, R. W.
and Gilliland, R. L., 2002, ApJ, 568, 377
Aumont,J.&Mac´ıas-P´erez,J.,F.,2007,MNRAS,
376, 739 Charbonneau,D.andAllen,L.E.,Megeath,S.T.,
Torres,G.,Alonso,R.,Brown,T.M.,etal.,A.,
Baluev, R. V., arXiv: 1209.3154v1
2005, ApJ, 626,523
10