# FPScausal: Functional propensity score weighting for causal inference with functional treatments, covariates, and outcomes

[![arXiv](https://img.shields.io/badge/arXiv-2608.03200-b31b1b.svg)](https://arxiv.org/abs/2608.03200)

`FPScausal` implements the *Functional Propensity Score (FPS) weighting*
methodology for causal inference with functional treatments and outcomes. If you use this package, please cite:

```
Ciardulli S. \& Fontana, N., Vantini S., Ieva, F. (2026). Generalized propensity score weighting for functional causal inference framework. arXiv. https://arxiv.org/abs/2608.03200. 
```

The package handles:

- **Functional treatment** 
- **Scalar and/or functional covariates**
- **Scalar, binary, or functional outcomes**

## Installation

Install the released version from CRAN:

```r
install.packages("FPScausal")
```

Or install the development version from GitHub:

```r
install.packages("devtools") # Install devtools if not already installed
devtools::install_github("NicoleFontana/FPSCausal")
```

## Methodology

1. **FPCA decomposition**: the functional treatment X(s) is represented via
   its Karhunen–Loève expansion, retaining the first L FPC scores.
2. **Empirical-likelihood balancing**: covariate-balancing weights are
   estimated by maximising the empirical likelihood subject to constraints
   that balance the FPC scores of the treatment against the observed
   confounders (and their interactions). The resulting dual problem is a
   smooth unconstrained minimisation solved via the BFGS quasi-Newton
   algorithm.
3. **Weighted least squares**: the causal effect function μ(s) (scalar
   outcome) or causal effect surface μ(s, t) (functional outcome) is
   recovered via weighted regression.
4. **Bootstrap CIs**: residual bootstrap (scalar/binary) or pairs bootstrap
   (functional outcome).

## Quick start

### Scalar outcome

```r
library(FPScausal)

# Simulate data — scalar covariates only
dat <- simulate_fps_data(
  n                    = 200,
  setting              = "LL",
  outcome_type         = "scalar",
  include_functional_cov = FALSE,
  seed                 = 42
)

# Step 1: estimate weights (treat_domain inferred from treat_grid)
w <- fps_weighting(
  treatment  = dat$X,
  treat_grid = dat$t_grid,
  covariates = dat$C
)

# Diagnostics
plot(w, type = "balance")
plot(w, type = "fpca_treatment")
plot(w, type = "weights")

# Step 2: estimate causal effect with bootstrap CIs
eff <- fps_effect_estimation(
  outcome    = dat$Y,
  fps_object = w,
  bootstrap  = TRUE,
  B          = 500,
  true_beta  = dat$true_beta,
  seed       = 1
)

plot(eff, type = "effect")       # μ(s) with CI ribbon and legend
plot(eff, type = "comparison")   # weighted vs unweighted
plot(eff, type = "significance") # significant time points
```

### Functional outcome

```r
dat_fn <- simulate_fps_data(
  n                    = 200,
  setting              = "LL",
  outcome_type         = "functional",
  include_functional_cov = FALSE,
  seed                 = 99
)

w_fn <- fps_weighting(
  treatment    = dat_fn$X,
  treat_grid   = dat_fn$t_grid,
  treat_domain = c(0, 1),
  domain_name  = "s",
  covariates   = dat_fn$C
)

plot(w_fn, type = "balance")
plot(w_fn, type = "fpca_treatment")

eff_fn <- fps_effect_estimation(
  outcome             = dat_fn$Y,
  fps_object          = w_fn,
  outcome_t_grid      = dat_fn$t_grid,
  outcome_domain      = c(0, 1),
  outcome_domain_name = "t",
  bootstrap           = TRUE,
  B                   = 500,
  seed                = 2
)

plot(eff_fn, type = "effect")           # μ(s,t) heatmap
plot(eff_fn, type = "bootstrap_slice",  # 1-D slice at t = 0.5
     point = 0.5, which_domain = "outcome")
plot(eff_fn, type = "bootstrap_slice",  # 1-D slice at s = 0.5
     point = 0.5, which_domain = "treatment")
plot(eff_fn, type = "significance")     # 2-D significance map
```

### Functional covariate 

```r
dat2 <- simulate_fps_data(
  n                    = 2000,
  setting              = "LL",
  outcome_type         = "scalar",
  include_functional_cov = TRUE,
  seed                 = 7
)

w2 <- fps_weighting(
  treatment   = dat2$X,
  treat_grid  = dat2$t_grid,
  domain_name = "s",
  covariates  = list(scalar = dat2$C, functional = list(dat2$D)),
  cov_grids   = list(dat2$t_grid)
)

plot(w2, type = "balance")
plot(w2, type = "fpca_covariates")

eff2 <- fps_effect_estimation(dat2$Y, w2, true_beta = dat2$true_beta)
plot(eff2, type = "effect")
```

## Package functions

| Function | Description |
|----------|-------------|
| `fps_weighting()` | Estimate FPS weights via empirical-likelihood balancing |
| `fps_effect_estimation()` | Estimate μ(s) or μ(s,t) with optional bootstrap CIs |
| `simulate_fps_data()` | Generate synthetic datasets (4 simulation settings) |
| `plot.fps_weighting()` | Balance, FPCA, and weight plots |
| `plot.fps_effect_estimation()` | Effect, comparison, slice, and significance plots |

## Simulation settings

`simulate_fps_data()` supports four settings varying whether the
treatment–confounder and confounder–outcome relationships are linear (L) or
nonlinear (N):

| Setting | Treatment–Confounder | Confounder–Outcome |
|---------|---------------------|--------------------|
| LL      | Linear              | Linear             |
| LN      | Linear              | Nonlinear          |
| NL      | Nonlinear           | Linear             |
| NN      | Nonlinear           | Nonlinear          |

## Dependencies

`fda`, `ggplot2`, `tidyr`, `MASS`, `wCorr`, `patchwork`, `progress`

## Reference

Ciardulli, S. and Fontana, N., Vantini S., Ieva F. (2026). *Functional propensity score weighting for causal inference with functional treatments, covariates, and outcomes*. arXiv:2608.03200. https://arxiv.org/abs/2608.03200

```bibtex
@misc{ciardulli2026generalizedpropensityscoreweighting,
      title={Generalized propensity score weighting for functional causal inference framework}, 
      author={Simone Ciardulli and Nicole Fontana and Simone Vantini and Francesca Ieva},
      year={2026},
      eprint={2608.03200},
      archivePrefix={arXiv},
      primaryClass={stat.ME},
      url={https://arxiv.org/abs/2608.03200}, 
}
```

## License

MIT
