joachim.py

import numpy as np
from sklearn.linear_model import LinearRegression

# Generate some example data
np.random.seed(42)
n_samples = 100

# Create synthetic features (x1: thickness, x2: stiffness, x3: scent, x4: another property)
X = np.random.rand(n_samples, 4) * 10  # Random values between 0 and 10

# Create synthetic target with known coefficients: y = x1^1 * x2^0 * x3^-1 * x4^0.5
true_coefficients = np.array([1.0, 0.0, -1.0, 0.5])
y = np.prod(X ** true_coefficients.reshape(1, -1), axis=1)

# Add some noise
y *= np.random.lognormal(0, 0.1, n_samples)

# Take log of all values (add small epsilon to avoid log(0))
epsilon = 1e-10
log_X = np.log(X + epsilon)
log_y = np.log(y + epsilon)

# Fit linear regression
reg = LinearRegression()
reg.fit(log_X, log_y)

# Get the coefficients and R-squared score
coefficients = reg.coef_
r2_score = reg.score(log_X, log_y)

print("True coefficients:", true_coefficients)
print("Fitted coefficients:", coefficients)
print("R-squared score:", r2_score)

# Test the model with a specific example
test_point = np.array([[2.0, 3.0, 1.5, 4.0]])
log_prediction = np.sum(coefficients * np.log(test_point))
prediction = np.exp(log_prediction)

print("\nTest prediction:")
print(f"Input values (x1, x2, x3, x4):", test_point[0])
print(f"Predicted y:", prediction)

# Calculate true value for comparison
true_value = np.prod(test_point ** true_coefficients)
print(f"True y:", true_value)