# Joint statistical modeling and deep learning for exoplanet detection and characterization by direct imaging at high contrast

Théo Bodrito, Olivier Flasseur, Julien Mairal, Jean Ponce, Maud Langlois, Anne-Marie Lagrange

► To cite this version:

Théo Bodrito, Olivier Flasseur, Julien Mairal, Jean Ponce, Maud Langlois, et al. Joint statistical modeling and deep learning for exoplanet detection and characterization by direct imaging at high contrast. EPSC-DPS Joint Meeting 2025, Sep 2025, Helsinki, Finland. pp.1-47. hal-05251737

HAL Id: hal-05251737

https://hal.science/hal-05251737v1

Submitted on 12 Sep 2025

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

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

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

Distributed under a Creative Commons Attribution 4.0 International License
# Joint statistical modeling and deep learning for exoplanet detection and characterization by direct imaging at high contrast

T. Bodrito¹, **O. Flasseur**², J. Mairal¹, J. Ponce¹, M. Langlois², A.-M. Lagrange³

multi-variate data

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

detection

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

characterization

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

reconstruction

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

<figure><img src="image_6.png" /><figcaption>Théo Bodrito</figcaption></figure>

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

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

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

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

<figure><img src="image_11.png" /></figure>
① Photons collection: large telescope (e.g., VLT)

② High angular resolution: adaptive optics (e.g., SPHERE)

③ High contrast ($\ge 10^6$): coronagraph (e.g., SPHERE)

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

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

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

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

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

BUT
exoplanets still
undetectable

④ Signal unmixing by dedicated algorithms
<figure><img src="image_17.png" /></figure>

PSF
=
planetary signal
<figure><img src="image_18.png" /></figure>

spatio-temporal slice cuts

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

PSF
=
planetary signal
<figure><img src="image_20.png" /></figure>

data

spatio-temporal slice cuts

PSF
=
planetary signal
<figure><img src="image_21.png" /></figure>

PSF
=
planetary signal
<figure><img src="image_22.png" /></figure>

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

PSF
=
planetary signal
Typical high-contrast data – VLT/SPHERE imager

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

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

PSF
=
planetary signal
<figure><img src="image_26.png" /></figure>

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

PSF
=
planetary signal
<figure><img src="image_28.png" /></figure>

* Faint signal from the exoplanets

* Non-stationary, strong, and spatially correlated nuisance component

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

PSF
=
planetary signal
<figure><img src="image_30.png" /></figure>

* Faint signal from the exoplanets

* Non-stationary, strong, and spatially correlated nuisance component

* Large fluctuations (stellar leakages)
<figure><img src="image_31.png" /></figure>

* Faint signal from the exoplanets

* Non-stationary, strong, and spatially correlated nuisance component

* Large fluctuations (stellar leakages)
<figure><img src="image_32.png" /></figure>

$\rightsquigarrow y \in \mathbb{R}^{H \times W \times T \times L} : H \times W \approx 10^6 \text{ pixels}, T \approx 100 \text{ frames}, L \in [2, 40] \text{ channels}$
<figure><img src="image_33.png" /></figure>

$\rightsquigarrow y \in \mathbb{R}^{H \times W \times T \times L} : H \times W \approx 10^6 \text{ pixels}, T \approx 100 \text{ frames}, L \in [2, 40] \text{ channels}$

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

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

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

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

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

<figure><img src="image_39.png" /></figure>
Extracting linear features from $\mathbf{y}^{(j)} \in \mathbb{R}^p$: $A\mathbf{y}^{(j)} \sim \mathcal{N}(\hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j)$, st. $A \in \mathbb{R}^{m \times p}$, $m \le p$

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

$\{A_j\}_{j \in \theta}$ (learned from database)

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

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

linear encoding

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

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

$\{\hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j\}_{j \in \chi(\mathbf{y})}$ (estimated for each obs.)
Extracting linear features from $\mathbf{y}^{(j)} \in \mathbb{R}^p$: $A\mathbf{y}^{(j)} \sim \mathcal{N}(\hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j)$, st. $A \in \mathbb{R}^{m \times p}$, $m \le p$

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

$\{A_j\}_{j \in \theta}$ (learned from database)

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

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

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

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

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

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

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

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

$\{\hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j\}_{j \in \chi(\mathbf{y})}$ (estimated for each obs.)
Extracting linear features from $\mathbf{y}^{(j)} \in \mathbb{R}^p$: $A\mathbf{y}^{(j)} \sim \mathcal{N}(\hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j)$, st. $A \in \mathbb{R}^{m \times p}$, $m \le p$

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

$\{A_j\}_{j \in \theta}$ (learned from database)

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

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

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

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

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

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

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

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

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

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

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

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

$\{\hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j\}_{j \in \chi}(\mathbf{y})$ (estimated for each obs.)
Extracting linear features from $\mathbf{y}^{(j)} \in \mathbb{R}^p$: $A\mathbf{y}^{(j)} \sim \mathcal{N}(\hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j)$, st. $A \in \mathbb{R}^{m \times p}$, $m \le p$

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

$\{A_j\}_{j \in \theta}$ (learned from database)

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

linear encoding

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

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

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

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

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

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

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

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

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

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

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

$\{\hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j\}_{j \in \chi(\mathbf{y})}$ (estimated for each obs.)
Integrating distributions from multiple spatial scales: patch-size $p \in \{8^2, 16^2, 32^2, 64^2\}$ pixels
linear encoding

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

$\{A_j\}_{j \in \theta}$ (learned from database) $\{\hat{m}_j, \hat{C}_j\}_{j \in \chi(y)}$ (estimated for each obs.)
Integrating distributions from multiple spatial scales: patch-size $p \in \{8^2, 16^2, 32^2, 64^2\}$ pixels
linear encoding

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

$\{A_j\}_{j \in \theta}$ (learned from database) $\{\hat{m}_j, \hat{C}_j\}_{j \in \chi(y)}$ (estimated for each pixel)
Integrating distributions from multiple spatial scales: patch-size $p \in \{8^2, 16^2, 32^2, 64^2\}$ pixels

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

$\{A_j\}_{j \in \theta}$ (learned from database)

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

linear encoding

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

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

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

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

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

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

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

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

$\{\hat{m}_j, \hat{C}_j\}_{j \in \chi(y)}$ (estimated for each obs.)
Accounting for rotational symmetries: $A[\mathbb{R}_{2\pi n/N}(\mathbf{y}_t)^{(j)}]_{n=0:N-1} \sim \mathcal{N}(\mathbf{m}_j, \mathbf{C}_j)$ st. $N = 1, 2, 4$

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

$\{A_j\}_{j \in \theta}$ (learned from database)

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

linear encoding

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

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

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

$\{\hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j\}_{j \in \chi(\mathbf{y})}$ (estimated for each obs.)
Accounting for rotational symmetries: $A[\mathbf{R}_{2\pi n/N}(\mathbf{y}_t)^{(j)}]_{n=0:N-1} \sim \mathcal{N}(\mathbf{m}_j, \mathbf{C}_j)$ st. $N = 1, 2, 4$

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

$\{A_j\}_{j \in \theta}$ (learned from database)

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

linear encoding

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

$A_j$

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

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

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

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

$A_j$

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

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

$\{\hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j\}_{j \in \chi(\mathbf{y})}$ (estimated for each obs.)
Accounting for rotational symmetries: $A[\mathbf{R}_{2\pi n/N}(\mathbf{y}_t)^{(j)}]_{n=0:N-1} \sim \mathcal{N}(\mathbf{m}_j, \mathbf{C}_j)$ st. $N = 1, 2, 4$

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

$\{A_j\}_{j \in \theta}$ (learned from database)

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

linear encoding

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

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

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

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

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

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

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

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

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

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

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

$\{\hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j\}_{j \in \chi(\mathbf{y})}$ (estimated for each obs.)
Accounting for rotational symmetries: $A[\mathbb{R}_{2\pi n/N}(\mathbf{y}_t)^{(j)}]_{n=0:N-1} \sim \mathcal{N}(\mathbf{m}_j, \mathbf{C}_j)$ st. $N=1,2,4$

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

$\{A_j\}_{j \in \theta}$ (learned from database)

$\mathbf{y}^{(j)}$

linear encoding

$A_j$

$\mathcal{N}(\hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j)$

$\{\hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j\}_{j \in \chi(\mathbf{y})}$ (estimated for each obs.)
Accounting for rotational symmetries: $A[\mathbb{R}_{2\pi n/N}(\mathbf{y}_t)^{(j)}]_{n=0:N-1} \sim \mathcal{N}(\mathbf{m}_j, \mathbf{C}_j)$ st. $N = 1, 2, 4$

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

$\{A_j\}_{j \in \theta}$ (learned from database)

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

$y^{(j)}$

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

$y^{(j)}$

linear encoding

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

$A_j$

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

$T$

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

$\mathcal{N}(\hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j)$

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

$A_j$

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

$T$

$\{\hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j\}_{j \in \chi(\mathbf{y})}$ (estimated for each obs.)
Accounting for rotational symmetries: $A[\mathbb{R}_{2\pi n/N}(\mathbf{y}_t)^{(j)}]_{n=0:N-1} \sim \mathcal{N}(\mathbf{m}_j, \mathbf{C}_j)$ st. $N = 1, 2, 4$

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

$\{A_j\}_{j \in \theta}$ (learned from database)

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

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

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

linear encoding

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

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

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

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

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

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

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

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

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

$\{\hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j\}_{j \in \chi(\mathbf{y})}$ (estimated for each obs.)
Data fidelity term: co-log-likelihood of the nuisance model

$$ \hat{\alpha} \in \underset{\alpha}{\operatorname{argmax}} \ell(\alpha, \mathbf{x}_0) := \prod_{t=1}^{T} \prod_{j \in S(\mathbf{x}_t)} \mathbb{P} \left( \mathbf{A}_j (\mathbf{y}_t^{(j)} - \alpha \mathbf{h}^{(j)}(\mathbf{x}_t)) \mid \hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j \right)^{w_j} $$

$S(\mathbf{x}_t)$: subset of modeled distributions around location $\mathbf{x}_t = r(\mathbf{x}_{t=1}, \phi_t)$ with learnable weight $w_j \in \mathbb{R}_+$ so that $\sum_j w_j = 1$

Learnable parameters: $\theta = \{\mathbf{A}_j, w_j\}_j$

Estimated parameters: $\chi = \{\hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j\}_j$ features' sample mean & shrunk covariance
**Data fidelity term: co-log-likelihood of the nuisance model**

$$ \hat{\alpha} \in \underset{\alpha}{\operatorname{argmax}} \ell(\alpha, \mathbf{x}_0) := \prod_{t=1}^{T} \prod_{j \in S(\mathbf{x}_t)} \mathbb{P} \left( \mathbf{A}_j (\mathbf{y}_t^{(j)} - \alpha \mathbf{h}^{(j)}(\mathbf{x}_t)) \middle| \hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j \right)^{w_j} $$

$S(\mathbf{x}_t)$: subset of modeled distributions around location $\mathbf{x}_t = r(\mathbf{x}_{t=1}, \phi_t)$ with learnable weight $w_j \in \mathbb{R}_+$ so that $\sum_j w_j = 1$

Learnable parameters: $\theta = \{\mathbf{A}_j, w_j\}_j$

Estimated parameters: $\chi = \{\hat{\mathbf{m}}_j, \hat{\mathbf{C}}_j\}_j$ features' sample mean & shrunk covariance

**Closed-form solution: flux estimator and detection criterion**

$$ \hat{\alpha}_{\lambda} = \frac{b_{\lambda}}{a_{\lambda}} \quad \text{and} \quad \hat{\sigma}_{\alpha_{\lambda}} = \frac{1}{\sqrt{a_{\lambda}}} \quad \Rightarrow \quad \gamma = \frac{\sum_{\lambda} b_{\lambda}}{\sqrt{\sum_{\lambda} a_{\lambda}}} \stackrel{\mathcal{H}_0}{\gtrsim} \tau, $$

$$ \text{with} = \begin{cases} b_{\lambda} = \sum_{t=1}^{T} \sum_{j \in S(\mathbf{x}_t)} w_j \mathbf{h}_{\lambda}^{(j)}(\mathbf{x}_t)^{\top} \mathbf{A}_j^{\top} \hat{\mathbf{C}}_j^{-1} (\mathbf{A}_j \mathbf{y}_{t,\lambda} - \hat{\mathbf{m}}_{\lambda,j}), & \lambda = 1, \dots, T \\ a_{\lambda} = \sum_{t=1}^{T} \sum_{j \in S(\mathbf{x}_t)} w_j \mathbf{h}_{\lambda}^{(j)}(\mathbf{x}_t)^{\top} \mathbf{A}_j^{\top} \hat{\mathbf{C}}_j^{-1} \mathbf{A}_j \mathbf{h}_{\lambda}^{(j)}(\mathbf{x}_t). & \lambda = 2, \dots, T \end{cases} $$
$$ egin{cases} b_{\lambda} = \sum_{t=1}^{T} \sum_{j \in S(\mathbf{x}_t)} w_j \mathbf{h}_{\lambda}^{(j)}(\mathbf{x}_t)^{\top} \mathbf{A}_j^{\top} \hat{\mathbf{C}}_j^{-1} (\mathbf{A}_j \mathbf{y}_{t,\lambda}^{(j)} - \hat{m}_{\lambda,j}), \\ a_{\lambda} = \sum_{t=1}^{T} \sum_{j \in S(\mathbf{x}_t)} w_j \mathbf{h}_{\lambda}^{(j)}(\mathbf{x}_t)^{\top} \mathbf{A}_j^{\top} \hat{\mathbf{C}}_j^{-1} \mathbf{A}_j \mathbf{h}_{\lambda}^{(j)}(\mathbf{x}_t). \end{cases} $$

⇒ learnable module integrated in a neural network architecture

<figure><img src="image_141.png" /></figure>
**Exoplanet-signal prior:**

* Learn a representation of exoplanet signatures

* Denoising flux estimation from the nuisance module with a U-Net

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

**Training:**

* Simulations of synthetic exoplanets in real data via direct model

* Full architecture ($\theta = \{\{A_j, w_j\}_{j,j}, \nu\}$) trained end-to-end

* Co-log-likelihood as loss function: $\mathcal{L}(\tilde{\alpha}, \tilde{\sigma}, \alpha_{\text{gt}}) = (\tilde{\alpha} - \alpha_{\text{gt}})/(2\tilde{\sigma}^2) + \log \tilde{\sigma}$
Detection performance: detection maps

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

min max △ false negative □ false positive min max △ false positive □ false negative □ false positive
Contrast curves: minimum detectable flux as a function of distance to the star

$$ \textit{Good observation} (\Delta_{\phi} = |\phi_{t=T} - \phi_{t=1}| = 75^{\circ}) $$

PACO (local statistical model)

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

proposed

<figure><img src="image_145.png" /></figure>
**Contrast curves:** minimum detectable flux as a function of distance to the star

*Medium observation* ($\Delta\phi = |\phi_{t=T} - \phi_{t=1}| = 25^{\circ}$)

PACO (local statistical model)

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

proposed

<figure><img src="image_147.png" /></figure>
**Contrast curves:** minimum detectable flux as a function of distance to the star

$$ \textit{Bad observation} (\Delta_{\phi} = |\phi_{t=T} - \phi_{t=1}| = 5^{\circ}) $$

PACO (local statistical model)

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

proposed

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

⇒ better detection sensitivity, allowing to process challenging datasets
* Hybrid modeling of strong & correlated nuisance component → joint multi-spectral processing with multiple scales and symmetries

* Unsupervised approach trained end-to-end

* Optimal fusion of multiple heterogeneous observations

* Reliable detection scores and estimates & improved sensitivity

## Perspectives:

* *Other modalities & instruments:* e.g. JWST with strong nuisance

* *Other tasks:* e.g. reconstruction of extended objects as circumstellar disks

<figure><img src="image_150.png" /></figure>
Thank you for your attention
Data acquisition

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

Post-processing

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

Patches extraction

Multi-scale

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

Rotational symmetries

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

Mixture of distributions

Convolutional statistical model

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

Detection by denoising

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

Learnable parameters

Estimated parameters

Hadamard product
Detection performance: ablation study

**Data:** 5 real observations from VLT/SPHERE instrument ($L = 2$)

**Metric:** precision-recall trade-off via AUC under ROCs (by varying $\tau$)

<table><thead><tr><th rowspan="2">Spatial scales</th><th colspan="3">Convolutional statistical model</th><th rowspan="2">+Learned denoiser</th></tr><tr><th>no symmetries</th><th>symmetries 180°</th><th>symmetries 180°, 90°</th></tr></thead><tbody><tr><td>8 × 8</td><td>0.554 ± 0.005</td><td>0.571 ± 0.005</td><td>0.575 ± 0.005</td><td rowspan="4">0.645 ± 0.005</td></tr><tr><td>+16 × 16</td><td>0.561 ± 0.005</td><td>0.573 ± 0.004</td><td>0.579 ± 0.005</td></tr><tr><td>+32 × 32</td><td>0.567 ± 0.005</td><td>0.575 ± 0.005</td><td>0.580 ± 0.005</td></tr><tr><td>+64 × 64</td><td>0.569 ± 0.005</td><td>0.577 ± 0.005</td><td>0.581 ± 0.005</td></tr><tr><td>8 × 8</td><td>0.713 ± 0.005</td><td>0.719 ± 0.005</td><td>0.723 ± 0.005</td><td rowspan="4">0.761 ± 0.004</td></tr><tr><td>+16 × 16</td><td>0.720 ± 0.005</td><td>0.723 ± 0.005</td><td>0.725 ± 0.005</td></tr><tr><td>+32 × 32</td><td>0.720 ± 0.005</td><td>0.725 ± 0.005</td><td>0.726 ± 0.005</td></tr><tr><td>+64 × 64</td><td>0.720 ± 0.005</td><td>0.725 ± 0.004</td><td>0.726 ± 0.004</td></tr></tbody></table>

⇒ Joint multi-spectral processing, multiple spatial scales, rotational symmetries, and learned denoiser are beneficial
**Data: VLT/SPHERE IRDIS**

<table><thead><tr><th>Modality</th><th>Method</th><th>HD159911 (54°)</th><th>HD216803 (23°)</th><th>HD206860 (11°)</th><th>HD188228 (6°)</th><th>HD102647 (2°)</th><th>average AUC</th></tr></thead><tbody><tr><td rowspan="5">ADI</td><td>cADI</td><td>0.288 ± 0.014</td><td>0.422 ± 0.007</td><td>0.489 ± 0.010</td><td>0.303 ± 0.014</td><td>0.343 ± 0.009</td><td>0.369 ± 0.005</td></tr><tr><td>PCA</td><td>0.634 ± 0.010</td><td>0.643 ± 0.011</td><td>0.505 ± 0.009</td><td>0.392 ± 0.011</td><td>0.218 ± 0.011</td><td>0.478 ± 0.005</td></tr><tr><td>PACO</td><td>0.629 ± 0.006</td><td>0.669 ± 0.012</td><td>0.579 ± 0.015</td><td>0.517 ± 0.015</td><td>0.207 ± 0.012</td><td>0.520 ± 0.006</td></tr><tr><td>MODEL&amp;CO</td><td>0.653 ± 0.010</td><td>0.731 ± 0.011</td><td>0.646 ± 0.012</td><td><strong>0.638</strong> ± 0.014</td><td><strong>0.554</strong> ± 0.009</td><td><strong>0.645</strong> ± 0.005</td></tr><tr><td>Proposed</td><td><strong>0.673</strong> ± 0.009</td><td><strong>0.740</strong> ± 0.010</td><td><strong>0.661</strong> ± 0.012</td><td>0.631 ± 0.013</td><td>0.518 ± 0.012</td><td><strong>0.645</strong> ± 0.005</td></tr><tr><td rowspan="4">ASDI</td><td>cASDI</td><td>0.448 ± 0.009</td><td>0.537 ± 0.013</td><td>0.451 ± 0.007</td><td>0.352 ± 0.017</td><td>0.294 ± 0.012</td><td>0.417 ± 0.006</td></tr><tr><td>PCA</td><td>0.694 ± 0.013</td><td>0.696 ± 0.007</td><td>0.552 ± 0.011</td><td>0.398 ± 0.011</td><td>0.236 ± 0.009</td><td>0.515 ± 0.005</td></tr><tr><td>PACO</td><td>0.698 ± 0.015</td><td>0.768 ± 0.008</td><td>0.710 ± 0.014</td><td>0.700 ± 0.009</td><td>0.589 ± 0.014</td><td>0.693 ± 0.005</td></tr><tr><td>Proposed</td><td><strong>0.731</strong> ± 0.009</td><td><strong>0.804</strong> ± 0.013</td><td><strong>0.747</strong> ± 0.010</td><td><strong>0.782</strong> ± 0.008</td><td><strong>0.744</strong> ± 0.006</td><td><strong>0.761</strong> ± 0.004</td></tr></tbody></table>

**Detection:** aggregated area under ROC curves (precision-recall trade-off)

<table><thead><tr><th>Method</th><th>flux error (ARE)</th><th>position error (RMSE)</th></tr></thead><tbody><tr><td>PACO</td><td>0.56</td><td>0.21</td></tr><tr><td>Proposed</td><td><strong>0.51</strong></td><td><strong>0.11</strong></td></tr></tbody></table>

**Astrometry & photometry:** absolute relative error (ARE) and root Mean Square Error (RMSE)
# Fusion of multiple observations of the same star

**Task:** joint detection and orbit ($\mu \in \mathbb{R}^7$) estimation from $N$ observations.

**Method:** MLE($\mu$) involves sufficient statistics $\{a_n, b_n\}_{n=1:N}$:

$$ \hat{\mu} = \underset{\mu}{\operatorname{argmax}} \left\{ \mathcal{C}(\mu) = \sum_{n=1}^{N} \frac{\left( \left[ b_n(x_n(\mu)) \right]_+ \right)^2}{a_n(x_n(\mu))} \right\} \Rightarrow \hat{\gamma}^{\text{comb}} \approx \sqrt{\mathcal{C}(\mu)}. $$

→ Solved via discrete exploration or using HMC algorithm (Dallant+ 2024)

**Example:** 23 observations
of the HR 8799 star
acquired between
2015 and 2022.

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

<figure><img src="image_158.png" /></figure>
<figure><img src="image_159.png" /></figure>

Results:

<figure><img src="image_160.png" /><figcaption>reconstructed speckles
(normalized by tiles)</figcaption></figure>

<figure><img src="image_161.png" /><figcaption>reconstructed object</figcaption></figure>

<figure><img src="image_162.png" /><figcaption>groundtruth object</figcaption></figure>

<figure><img src="image_163.png" /></figure>
Results: inferred detection maps – *parallactic rotation* ≈ 23°

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

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

<figure><img src="image_166.png" /></figure>
Results: inferred detection maps – parallactic rotation ≈ 10°

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

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

<figure><img src="image_169.png" /></figure>
<figure><img src="image_170.png" /><figcaption>Results: inferred detection maps – *parallactic rotation* ≈ 6°</figcaption></figure>

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

<figure><img src="image_172.png" /></figure>
<figure><img src="image_173.png" /><figcaption>Results: inferred detection maps – *parallactic rotation* ≈ 2°!</figcaption></figure>

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

<figure><img src="image_175.png" /></figure>
* Low nb of samples ($T \approx L \approx K$) ⇒ $\hat{S}_n^{\text{spat}}$ and $\hat{S}_n^{\text{spec}}$ noisy/rank-deficient

* Data-driven and spatially adaptive shrinkage:

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

* Optimal setting through risk minimization (Oracle Shrinkage):

$$ \hat{\rho}^{\text{OS}} = \arg \min_{\rho \in [0,1]} \frac{\mathbb{E}(\|\hat{\mathbf{C}} - \mathbf{C}\|_{\text{F}}^2)}{\gamma \mathbb{E}(\text{tr}((\hat{\mathbf{F}} - \hat{\mathbf{S}})^2))}, \quad \text{with } \mathbb{E}(\hat{\mathbf{S}}) = \gamma^{-1} \mathbf{C}. $$

* Practical closed-form (Oracle Approximate Shrinkage):

$$ \hat{\rho}^{\text{OAS}} = \frac{(\gamma\nu + \gamma^{-1})(\text{tr}(\hat{\mathbf{S}}^2) - \sum_i [\hat{\mathbf{S}}]_{ii}^2) + \gamma\eta(\text{tr}^2(\hat{\mathbf{S}}) - \sum_i [\hat{\mathbf{S}}]_{ii}^2)}{\gamma\nu(\text{tr}(\hat{\mathbf{S}}^2) - \sum_i [\hat{\mathbf{S}}]_{ii}^2)}. $$