---
title: "Use Case: Genomic Data Smoothing"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Use Case: Genomic Data Smoothing}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

<!-- markdownlint-disable MD024 -->

```{r setup, include=FALSE}
knitr::opts_chunk$set(collapse = TRUE, comment = "#>")
```

## Overview

LOESS is well-suited for genomic data such as DNA methylation profiles, ChIP-seq
signals, and expression arrays.

---

## Methylation Profile Smoothing

### The Challenge

DNA methylation data (from bisulfite sequencing or arrays) shows
position-dependent patterns that can be obscured by measurement noise.

### Solution

A small `fraction = 0.1` lets LOESS follow fine-scale spatial structure without
smearing the transitions between methylated and unmethylated regions.
`intervals = intervals_opts(confidence = 0.95)` produces uncertainty bands that naturally widen at
positions with sparser CpG coverage, making low-confidence segments immediately
apparent in the plot.

```{r use_case_genomics_1}
library(rfastloess)

set.seed(42)
n <- 1000
positions <- sort(runif(n, 0, 1e6))

true_meth <- 0.5 + 0.3 * sin(positions / 1e5)
observed  <- true_meth + rnorm(n, sd = 0.15)
observed  <- pmax(0, pmin(1, observed))

model <- Loess(
    fraction = 0.1,
    iterations = 3,
    intervals = intervals_opts(confidence = 0.95)
)
result <- fit(model, positions, observed)

plot(positions, observed, pch = ".", col = "gray",
    xlab = "Genomic Position (bp)", ylab = "Methylation Level",
    main = "Methylation Profile Smoothing")
lines(result$x, result$y, col = "blue", lwd = 2)
lines(result$x, result$confidence_lower, col = "blue", lty = 2)
lines(result$x, result$confidence_upper, col = "blue", lty = 2)
legend("topright", c("Observed", "Smoothed", "95% CI"),
        pch = c(1, NA, NA), lty = c(NA, 1, 2),
        col = c("gray", "blue", "blue"))
```

---

## ChIP-seq Signal Smoothing

### Application

ChIP-seq experiments produce sparse, noisy coverage data. LOESS can help
identify binding regions.

`fraction = 0.05` provides high spatial resolution — important for resolving
narrow binding peaks that would otherwise be smeared into the background. The
larger `iterations = 5` is deliberate: Poisson-distributed read counts
produce tall, isolated spikes, and extra robustness iterations progressively
down-weight them so the estimated background level is not inflated by a
handful of extreme counts.

```{r use_case_genomics_2}
library(rfastloess)
set.seed(123)

positions <- seq(0, 9990, by = 10)
n <- length(positions)

background <- 10
peak1 <- 50 * exp(-((positions - 2000)^2) / (2 * 200^2))
peak2 <- 80 * exp(-((positions - 5000)^2) / (2 * 300^2))
peak3 <- 40 * exp(-((positions - 8000)^2) / (2 * 150^2))

true_signal <- background + peak1 + peak2 + peak3
observed <- rpois(n, true_signal)

model <- Loess(
    fraction = 0.05,   # Very local smoothing,
    iterations = 5,    # Strong robustness,
    outputs = "residuals"
)
result <- fit(model, positions, observed)

# Identify peaks (smoothed signal significantly above background)
threshold <- quantile(result$y, 0.75)
peaks <- positions[result$y > threshold]
cat(sprintf("Found %d peak positions (first 5): %s\n",
            length(peaks), toString(head(peaks, 5))))
```

---

## Large Genome Coverage (Streaming)

For whole-genome data that doesn't fit in memory:

```{r use_case_genomics_3}
library(rfastloess)
set.seed(42)

positions <- seq(0, 9990, by = 10)
coverage  <- rpois(length(positions), 50)

# Process chromosome-by-chromosome or in chunks
model <- StreamingLoess(
    fraction   = 0.05,
    chunk_size = 100000,   # 100kb chunks,
    overlap    = 10000,    # 10kb overlap,
    merge_strategy = "weighted_average"
)
process_chunk(model, positions, coverage)
result <- finalize(model)
cat("Smoothed y[0]:", result$y[1], "\n")
```

---

## Best Practices for Genomic Data

| Consideration | Recommendation |
| --- | --- |
| **Fraction** | 0.05–0.15 (preserve local features) |
| **Iterations** | 3–5 (handle sequencing outliers) |
| **Large data** | Use streaming mode |
| **Sparse regions** | Use `boundary_policy = "extend"` |
| **Multiple chromosomes** | Process separately or ensure sorted |

---

## See Also

- `vignette("concepts", package = "rfastloess")` — How LOESS works
- `?Loess` — All options
- `vignette("robustness", package = "rfastloess")` — Outlier
  downweighting in depth
- `vignette("merge", package = "rfastloess")` — Streaming chunk reconciliation
- `vignette("boundary", package = "rfastloess")` — Edge handling for
  sparse regions
- `vignette("use-case-real-time", package = "rfastloess")` — For sequencing runs

```{r session-info}
sessionInfo()
```