Skip to content

Commit 3c9cb1b

Browse files
committed
Multiple polynomial support and nan checking added
1 parent a9df5b1 commit 3c9cb1b

1 file changed

Lines changed: 66 additions & 44 deletions

File tree

codonPython/tolerance.py

Lines changed: 66 additions & 44 deletions
Original file line numberDiff line numberDiff line change
@@ -6,7 +6,7 @@
66
from statsmodels.sandbox.regression.predstd import wls_prediction_std
77

88

9-
def check_tolerance(t, y, to_exclude: int = 1, poly_features: int = 2, alpha: float = 0.05) -> pd.DataFrame:
9+
def check_tolerance(t, y, to_exclude: int = 1, poly_features: list = [1, 2], alpha: float = 0.05) -> pd.DataFrame:
1010
"""
1111
Check that some future values are within a weighted least squares confidence interval.
1212
@@ -18,8 +18,10 @@ def check_tolerance(t, y, to_exclude: int = 1, poly_features: int = 2, alpha: fl
1818
The corresponding response variable values to X, of shape (N, 1).
1919
to_exclude : int, default = 1
2020
How many of the last y values will have their tolerances checked.
21-
poly_features : int, default = 2
22-
Degree of polynomial features to fit to the data.
21+
poly_features : list, default = [1, 2]
22+
List of degrees of polynomial features to fit to the data. One model will be
23+
produced for each number in the list, eg. the default will fit a linear and
24+
a second degree polynomial to the data and return both sets of results.
2325
alpha : float, default = 0.05
2426
Alpha parameter for the weighted least squares confidence interval.
2527
@@ -28,61 +30,81 @@ def check_tolerance(t, y, to_exclude: int = 1, poly_features: int = 2, alpha: fl
2830
-------
2931
pd.DataFrame
3032
DataFrame of shape (to_exclude, 4) containing:
31-
"yhat_u" : Upper condfidence interval for y
32-
"yobs" : Observed value for y
33-
"yhat" : Predicted value for y
34-
"yhat_l" : Lower confidence interval for y
33+
"yhat_u" : Upper condfidence interval for y
34+
"yobs" : Observed value for y
35+
"yhat" : Predicted value for y
36+
"yhat_l" : Lower confidence interval for y
37+
"polynomial": Degree of polynomial fit to data
3538
3639
3740
Examples
3841
--------
3942
>>> check_tolerance(
4043
... t = np.array([1001,1002,1003,1004,1005,1006]),
4144
... y = np.array([2,3,4,4.5,5,5.1]),
42-
... ).round(3).to_dict()
43-
{'yhat_u': {0: 6.061}, 'yobs': {0: 5.1}, 'yhat': {0: 5.2}, 'yhat_l': {0: 4.339}}
45+
... to_exclude = 2,
46+
... poly_features = [2],
47+
... )
48+
yhat_u yobs yhat yhat_l polynomial
49+
0 9.077182 5.0 4.875 0.672818 2
50+
1 13.252339 5.1 4.975 -3.302339 2
4451
"""
4552

46-
if not isinstance(poly_features, int) or 0 >= poly_features >= 4:
47-
raise ValueError("Please input an integer from 0 to 4 for poly_features.")
53+
if not isinstance(poly_features, list):
54+
raise ValueError("Please input a list of integers from 0 to 4 for poly_features.")
55+
assert all([0 <= degree <= 4 for degree in poly_features]), (
56+
"Please ensure all numbers in poly_features are from 0 to 4."
57+
)
4858
if not isinstance(alpha, float) or 0 >= alpha >= 1:
4959
raise ValueError("Please input a float between 0 and 1 for alpha.")
60+
if not isinstance(to_exclude, int) or len(t) <= to_exclude < 1:
61+
raise ValueError("Please input an integer between 1 and your sample size for to_exclude.")
62+
assert ((len(t) - to_exclude) >= 4), (
63+
"""The sample size for your model is smaller than 4. This will not produce a good
64+
model. Either reduce to_exclude or increase your sample size to continue."""
65+
)
66+
assert np.isfinite(y).all(), (
67+
"Your sample contains missing or infinite values for y. Please exclude these values to continue."
68+
)
5069

51-
N = len(t)
52-
53-
if not isinstance(to_exclude, int) or N <= to_exclude < 1:
54-
raise ValueError("Please input an integer between 1 and your sample size to exclude.")
55-
if N < 4:
56-
raise ValueError("Your sample size is smaller than 4. This will not produce a good model.")
5770

5871
# Sort data by X increasing
5972
idx = np.argsort(t)
6073
t = t[idx]
6174
y = y[idx]
62-
63-
transforms = make_pipeline(
64-
StandardScaler(),
65-
PolynomialFeatures(degree=poly_features),
66-
)
67-
68-
# Fit transforms to train data, apply them to all data
69-
fitted_transforms = transforms.fit(t[:-to_exclude].reshape(-1, 1))
70-
t = fitted_transforms.transform(t.reshape(-1, 1))
71-
72-
t_train, y_train = t[:-to_exclude, :], y[:-to_exclude]
73-
t_predict, y_predict = t[-to_exclude:, :], y[-to_exclude:]
74-
75-
# Fit ordinary least squares model to the training data, then predict for the
76-
# prediction data.
77-
model = sm.OLS(y_train, t_train).fit()
78-
yhat = model.predict(t_predict)
79-
80-
# Calculate confidence interval of fitted model.
81-
_, yhat_l, yhat_u = wls_prediction_std(model, t_predict, alpha=alpha)
82-
83-
return pd.DataFrame({
84-
"yhat_u" : yhat_u,
85-
"yobs" : y_predict,
86-
"yhat" : yhat,
87-
"yhat_l" : yhat_l,
88-
})
75+
results = pd.DataFrame()
76+
77+
for degree in poly_features:
78+
transforms = make_pipeline(
79+
StandardScaler(),
80+
PolynomialFeatures(degree=degree),
81+
)
82+
83+
# Fit transforms to training data, apply to all data.
84+
fitted_transforms = transforms.fit(t[:-to_exclude].reshape(-1, 1))
85+
_t = fitted_transforms.transform(t.reshape(-1, 1))
86+
87+
t_train, y_train = _t[:-to_exclude, :], y[:-to_exclude]
88+
t_predict, y_predict = _t[-to_exclude:, :], y[-to_exclude:]
89+
90+
# Fit ordinary least squares model to the training data, then predict for the
91+
# prediction data.
92+
model = sm.OLS(y_train, t_train).fit()
93+
yhat = model.predict(t_predict)
94+
95+
# Calculate confidence interval of fitted model.
96+
_, yhat_l, yhat_u = wls_prediction_std(model, t_predict, alpha=alpha)
97+
98+
# Store model results in master frame
99+
results = results.append(
100+
pd.DataFrame({
101+
"yhat_u" : yhat_u,
102+
"yobs" : y_predict,
103+
"yhat" : yhat,
104+
"yhat_l" : yhat_l,
105+
"polynomial" : degree
106+
})
107+
)
108+
109+
return results
110+

0 commit comments

Comments
 (0)