///|
/// Residual diagnostics for validating fitted robust models.
pub struct ResidualDiagnostic {
  name : String
  value : Double
  threshold : Double
  passed : Bool
  interpretation : String
}

///|
pub struct ResidualReport {
  count : Int
  center : Double
  scale : Double
  diagnostics : Array[ResidualDiagnostic]
  score : Double
}

///|
pub fn residual_diagnostic(
  name : String,
  value : Double,
  threshold : Double,
  passed : Bool,
  interpretation : String,
) -> ResidualDiagnostic {
  { name, value, threshold, passed, interpretation }
}

///|
pub fn residual_report(residuals : Array[Double]) -> ResidualReport {
  let diagnostics = []
  let center = mean(residuals)
  let scale = mad(residuals)
  let autocorrelation = abs_double(signal_autocorrelation(residuals, 1))
  let outlier_rate = if residuals.length() == 0 {
    0.0
  } else {
    outlier_indices_z(residuals, threshold=3.5).length().to_double() /
    residuals.length().to_double()
  }
  diagnostics.push(
    residual_diagnostic(
      "bias",
      abs_double(center),
      0.1,
      abs_double(center) <= 0.1,
      "residual center",
    ),
  )
  diagnostics.push(
    residual_diagnostic(
      "autocorrelation",
      autocorrelation,
      0.2,
      autocorrelation <= 0.2,
      "lag-one dependence",
    ),
  )
  diagnostics.push(
    residual_diagnostic(
      "outlier_rate",
      outlier_rate,
      0.05,
      outlier_rate <= 0.05,
      "extreme residual fraction",
    ),
  )
  diagnostics.push(
    residual_diagnostic(
      "scale",
      scale,
      1.0e12,
      scale <= 1.0e12,
      "residual dispersion",
    ),
  )
  let mut passed = 0
  for diagnostic in diagnostics {
    if diagnostic.passed {
      passed += 1
    }
  }
  {
    count: residuals.length(),
    center,
    scale,
    diagnostics,
    score: if diagnostics.length() == 0 {
      1.0
    } else {
      passed.to_double() / diagnostics.length().to_double()
    },
  }
}

///|
pub fn residual_diagnostics(
  data : Array[Double],
  fitted : Array[Double],
) -> ResidualReport {
  let residuals = []
  let count = if data.length() < fitted.length() {
    data.length()
  } else {
    fitted.length()
  }
  for index = 0; index < count; index = index + 1 {
    residuals.push(data[index] - fitted[index])
  }
  residual_report(residuals)
}

///|
pub fn residual_values(report : ResidualReport) -> Array[Double] {
  [report.count.to_double(), report.center, report.scale, report.score]
}

///|
pub fn residual_diagnostic_names(report : ResidualReport) -> Array[String] {
  let result = []
  for diagnostic in report.diagnostics {
    result.push(diagnostic.name)
  }
  result
}

///|
pub fn residual_diagnostic_scores(report : ResidualReport) -> Array[Double] {
  let result = []
  for diagnostic in report.diagnostics {
    result.push(diagnostic.value)
  }
  result
}

///|
pub fn residual_failed(report : ResidualReport) -> Array[ResidualDiagnostic] {
  let result = []
  for diagnostic in report.diagnostics {
    if !diagnostic.passed {
      result.push(diagnostic)
    }
  }
  result
}

///|
pub fn residual_report_lines(report : ResidualReport) -> Array[String] {
  let lines = [
    "count=" + report.count.to_string(),
    "center=" + report.center.to_string(),
    "scale=" + report.scale.to_string(),
    "score=" + report.score.to_string(),
  ]
  for diagnostic in report.diagnostics {
    lines.push(
      diagnostic.name +
      "=" +
      diagnostic.value.to_string() +
      " threshold=" +
      diagnostic.threshold.to_string() +
      " passed=" +
      diagnostic.passed.to_string(),
    )
  }
  lines
}

///|
pub fn residual_report_string(report : ResidualReport) -> String {
  residual_report_lines(report).join("\n")
}

///|
pub fn residual_is_acceptable(
  report : ResidualReport,
  threshold : Double,
) -> Bool {
  report.score >= threshold
}

///|
pub fn residual_normalized(data : Array[Double]) -> Array[Double] {
  transform_robust_scale(data)
}

///|
pub fn residual_studentized(data : Array[Double]) -> Array[Double] {
  let center = mean(data)
  let scale = sample_stddev(data)
  let result = []
  for value in data {
    result.push(if scale == 0.0 { 0.0 } else { (value - center) / scale })
  }
  result
}

///|
pub fn residual_rolling_bias(
  data : Array[Double],
  window : Int,
) -> Array[Double] {
  let result = []
  let width = if window < 1 { 1 } else { window }
  for index = 0; index < data.length(); index = index + 1 {
    let start = if index + 1 > width { index + 1 - width } else { 0 }
    let values = []
    for cursor = start; cursor <= index; cursor = cursor + 1 {
      values.push(data[cursor])
    }
    result.push(mean(values))
  }
  result
}

///|
pub fn residual_rolling_scale(
  data : Array[Double],
  window : Int,
) -> Array[Double] {
  let result = []
  let width = if window < 1 { 1 } else { window }
  for index = 0; index < data.length(); index = index + 1 {
    let start = if index + 1 > width { index + 1 - width } else { 0 }
    let values = []
    for cursor = start; cursor <= index; cursor = cursor + 1 {
      values.push(data[cursor])
    }
    result.push(mad(values))
  }
  result
}

///|
pub fn residual_whiten(data : Array[Double]) -> Array[Double] {
  let result = []
  if data.length() == 0 {
    return result
  }
  result.push(data[0] - mean(data))
  for index = 1; index < data.length(); index = index + 1 {
    result.push(data[index] - data[index - 1])
  }
  result
}

///|
pub fn residual_ljung_like(data : Array[Double], maximum_lag : Int) -> Double {
  let path = signal_autocorrelation_path(data, maximum_lag)
  let mut total = 0.0
  for index = 1; index < path.length(); index = index + 1 {
    total += path[index] * path[index]
  }
  total * data.length().to_double()
}

///|
pub fn residual_quantiles(data : Array[Double]) -> Array[Double] {
  [
    quantile(data, 0.01),
    quantile(data, 0.05),
    quantile(data, 0.25),
    quantile(data, 0.5),
    quantile(data, 0.75),
    quantile(data, 0.95),
    quantile(data, 0.99),
  ]
}

///|
pub fn residual_symmetry(data : Array[Double]) -> Double {
  let left = abs_double(quantile(data, 0.25) - median(data))
  let right = abs_double(quantile(data, 0.75) - median(data))
  if left + right == 0.0 {
    1.0
  } else {
    1.0 - abs_double(left - right) / (left + right)
  }
}

///|
pub fn residual_robustness(data : Array[Double]) -> Double {
  let normalized = residual_normalized(data)
  1.0 / (1.0 + mad(normalized))
}

///|
pub fn residual_compare(
  raw : Array[Double],
  cleaned : Array[Double],
) -> Array[Double] {
  let raw_report = residual_report(raw)
  let clean_report = residual_report(cleaned)
  [
    raw_report.score,
    clean_report.score,
    clean_report.score - raw_report.score,
    raw_report.scale,
    clean_report.scale,
  ]
}

///|
pub fn residual_summary(data : Array[Double]) -> Array[Double] {
  let report = residual_report(data)
  [
    report.count.to_double(),
    report.center,
    report.scale,
    report.score,
    signal_autocorrelation(data, 1),
    residual_symmetry(data),
    residual_robustness(data),
  ]
}