///|
pub struct DoseResponsePoint {
  dose : Double
  count : Int
  mean_outcome : Double
  standard_error : Double
}

///|
pub fn dose_response_curve(
  dose : Array[Double],
  outcomes : Array[Double],
  bins : Int,
) -> Array[DoseResponsePoint] {
  let count_bins = if bins < 1 { 1 } else { bins }
  let result : Array[DoseResponsePoint] = Array::new(capacity=count_bins)
  let minimum = quantile(dose, 0.0)
  let maximum = quantile(dose, 1.0)
  let width = if maximum == minimum {
    1.0
  } else {
    (maximum - minimum) / count_bins.to_double()
  }
  for bin in 0..= lower && dose[i] < upper {
        values.push(outcomes[i])
        doses.push(dose[i])
      }
    }
    result.push({
      dose: mean(doses),
      count: values.length(),
      mean_outcome: mean(values),
      standard_error: if values.length() <= 1 {
        0.0
      } else {
        std_dev(values) / values.length().to_double().sqrt()
      },
    })
  }
  result
}

///|
pub fn generalized_propensity_scores(
  dose : Array[Double],
  covariates : Array[Array[Double]],
) -> Array[Double] {
  let model = fit_linear_outcome_model(covariates, dose, ridge=1.0e-6)
  let predictions = predict_outcomes(model, covariates)
  let residuals = Array::new(capacity=dose.length())
  for i in 0.. Array[Double] {
  let n = if treatment.length() < propensity_scores.length() {
    treatment.length()
  } else {
    propensity_scores.length()
  }
  let result = Array::new(capacity=n)
  for i in 0.. Double {
  weighted_mean(outcomes, weights) + weighted_mean(dose, weights) * 0.0
}