Untitled

import numpy as np
import pandas as pd
from scipy.stats import linregress

# Assuming you have your data in a DataFrame named 'sensor_readings_t1'

# Define the list to store dictionaries for each linear regression model
linReg_CI_dicts = []

def sqrt_term(x_new, X_sensor, n):
    X_line = np.mean(X_sensor)
    S_xx = np.sum((X_sensor - X_line) ** 2) / (n - 1)

    if S_xx == 0:
        return 0.0

    return np.sqrt(1 / n + (x_new - X_line) ** 2 / S_xx)

def upper(x_new, a, b, Sigma_linReg, X_sensor, n, df):
    alpha = 0.95
    quantile = t.ppf(alpha, df)
    return x_new * a + b + quantile * Sigma_linReg * sqrt_term(x_new, X_sensor, n)

def lower(x_new, a, b, Sigma_linReg, X_sensor, n, df):
    alpha = 0.05
    quantile = t.ppf(alpha, df)
    return x_new * a + b - quantile * Sigma_linReg * sqrt_term(x_new, X_sensor, n)

def middle(x_new, a, b):
    return x_new * a + b




# Loop through each unique sensor
for x in range(7):
    # Select data for the current sensor
    data_sensor = sensor_readings_t1[sensor_readings_t1['Sensor'] == x]

    # Extract the dependent variable and independent variables for the current sensor
    y_sensor = np.sqrt(np.array(data_sensor['Estimated_VWC']))
    X_sensor = np.array(data_sensor['V_from_SWC'])


    # Add a constant to the independent variables matrix for the intercept term
    model = LinearRegression().fit(X_sensor.reshape(-1,1),y_sensor)

    # Get the confidence intervals for the coefficient estimates

    y_estimates = model.predict(X_sensor.reshape(-1,1))
    n = len(X_sensor)
    Sigma_linReg = np.sum((y_sensor - y_estimates)**2)/(n-2)

    # Degrees of freedom
    df = n-2


    a = copy.deepcopy(model.coef_[0])
    b = copy.deepcopy(model.intercept_)

    # Create the dictionary for the current sensor and add it to the list
    sensor_dict = {
        "Lower": lambda x, a=a, b=b, Sigma_linReg=Sigma_linReg, X_sensor=X_sensor, n=n, df=df: lower(x, a, b, Sigma_linReg, X_sensor, n, df),
        "Middle": lambda x, a=a,b=b: middle(x,a,b),
        "Upper": lambda x, a=a, b=b, Sigma_linReg=Sigma_linReg, X_sensor=X_sensor, n=n, df=df: upper(x, a, b, Sigma_linReg, X_sensor, n, df)
    }
    linReg_CI_dicts.append(sensor_dict)

# Now 'linReg_CI_dicts' contains dictionaries for each linear regression model with the functions to calculate the lower, middle, and upper bounds of the 95% confidence intervals.