Skip to content

Commit 771dd7f

Browse files
authored
Merge pull request #18 from codonlibrary/tolerance-checking
Tolerance checking added
2 parents 338c932 + f8a50db commit 771dd7f

4 files changed

Lines changed: 267 additions & 2 deletions

File tree

.travis.yml

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,6 @@
11
language: python
22
python:
3-
- 3.5
3+
- 3.6
44
install:
55
- pip install -r requirements.txt
66
- pip install codecov
Lines changed: 141 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,141 @@
1+
from codonPython.tolerance import check_tolerance
2+
import numpy as np
3+
import pandas as pd
4+
import pandas.util.testing as pdt
5+
import pytest
6+
7+
8+
testdata = [
9+
np.array([1234, 1235, 1236, 1237, 1238, 1239, 1240, 1241, 1242]),
10+
np.array([1, 2, 3, 4, 5, 5.5, 6, 6.5, 7]),
11+
]
12+
13+
14+
@pytest.mark.parametrize("t, y, to_exclude, poly_features, alpha, expected", [
15+
(
16+
*testdata,
17+
2,
18+
[1, 2],
19+
0.05,
20+
pd.DataFrame({
21+
'yhat_u': [
22+
8.11380197739608,
23+
9.051653693670929,
24+
7.127135023632205,
25+
7.735627110021585,
26+
],
27+
'yobs': [6.5, 7.0, 6.5, 7.0],
28+
'yhat': [
29+
7.214285714285714,
30+
8.071428571428573,
31+
6.500000000000002,
32+
6.821428571428574,
33+
],
34+
'yhat_l': [
35+
6.31476945117535,
36+
7.091203449186216,
37+
5.872864976367799,
38+
5.907230032835563,
39+
],
40+
'polynomial': [1, 1, 2, 2]
41+
})
42+
),
43+
(
44+
*testdata,
45+
2,
46+
[3],
47+
0.05,
48+
pd.DataFrame({
49+
'yhat_u': [
50+
6.753927165005773,
51+
7.214574732953706,
52+
],
53+
'yobs': [6.5, 7.0],
54+
'yhat': [
55+
6.0000000000000036,
56+
5.571428571428576,
57+
],
58+
'yhat_l': [
59+
5.2460728349942345,
60+
3.928282409903445,
61+
],
62+
'polynomial': [3, 3]
63+
})
64+
),
65+
])
66+
def test_tolerance_checking_BAU(t, y, to_exclude, poly_features, alpha, expected):
67+
obtained = check_tolerance(
68+
t,
69+
y,
70+
to_exclude=to_exclude,
71+
poly_features=poly_features,
72+
alpha=alpha,
73+
)
74+
pdt.assert_frame_equal(expected, obtained)
75+
76+
77+
@pytest.mark.parametrize("t, y, to_exclude, poly_features, alpha", [
78+
(
79+
*testdata,
80+
2,
81+
"flamingo", # This should be a list
82+
0.05,
83+
),
84+
(
85+
*testdata,
86+
2,
87+
[2],
88+
"flamingo", # Needs to be int
89+
),
90+
(
91+
*testdata,
92+
2,
93+
[2],
94+
42, # Needs to be between 0 and 1
95+
),
96+
(
97+
*testdata,
98+
"flamingo", # Needs to be int
99+
[2],
100+
0.05,
101+
),
102+
])
103+
def test_ValueErrors(t, y, to_exclude, poly_features, alpha):
104+
with pytest.raises(ValueError):
105+
check_tolerance(t, y, to_exclude=to_exclude,
106+
poly_features=poly_features, alpha=alpha)
107+
108+
109+
@pytest.mark.parametrize("t, y, to_exclude, poly_features, alpha", [
110+
(
111+
*testdata,
112+
2,
113+
[42], # Elements in the list should be between 0 and 4
114+
0.05,
115+
),
116+
(
117+
*testdata,
118+
42, # Can't have to_exclude making your sample size smaller than 4
119+
[2],
120+
0.05,
121+
),
122+
(
123+
np.array([1234, 1235, 1236, 1237, 1238, 1239,
124+
1240, 1241, np.nan]), # Missing t value
125+
np.array([1, 2, 3, 4, 5, 5.5, 6, 6.5, 7]),
126+
2,
127+
[2],
128+
0.05,
129+
),
130+
(
131+
np.array([1234, 1235, 1236, 1237, 1238, 1239, 1240, 1241, 1242]),
132+
np.array([1, 2, 3, 4, 5, 5.5, 6, 6.5, np.nan]), # Missing y value
133+
2,
134+
[2],
135+
0.05,
136+
)
137+
])
138+
def test_AssertionErrors(t, y, to_exclude, poly_features, alpha):
139+
with pytest.raises(AssertionError):
140+
check_tolerance(t, y, to_exclude=to_exclude,
141+
poly_features=poly_features, alpha=alpha)

codonPython/tolerance.py

Lines changed: 122 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,122 @@
1+
import numpy as np
2+
import pandas as pd
3+
from sklearn.preprocessing import StandardScaler, PolynomialFeatures
4+
from sklearn.pipeline import make_pipeline
5+
import statsmodels.api as sm
6+
from statsmodels.sandbox.regression.predstd import wls_prediction_std
7+
8+
9+
def check_tolerance(t, y, to_exclude: int = 1, poly_features: list = [1, 2], alpha: float = 0.05, predict_all: bool = False) -> pd.DataFrame:
10+
"""
11+
Check that some future values are within a weighted least squares confidence interval.
12+
13+
Parameters
14+
----------
15+
t : np.array
16+
N explanatory time points of shape (N, 1).
17+
y : np.array
18+
The corresponding response variable values to X, of shape (N, 1).
19+
to_exclude : int, default = 1
20+
How many of the last y values will have their tolerances checked.
21+
poly_features : list, default = [1, 2]
22+
List of degrees of polynomial basis 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.
25+
alpha : float, default = 0.05
26+
Alpha parameter for the weighted least squares confidence interval.
27+
predict_all : bool, default = False
28+
Set to true to show predictions for all points of the dataset.
29+
30+
31+
Returns
32+
-------
33+
pd.DataFrame
34+
DataFrame containing:
35+
"yhat_u" : Upper condfidence interval for y
36+
"yobs" : Observed value for y
37+
"yhat" : Predicted value for y
38+
"yhat_l" : Lower confidence interval for y
39+
"polynomial": Max polynomial of model fit to the data
40+
41+
42+
Examples
43+
--------
44+
>>> check_tolerance(
45+
... t = np.array([1001,1002,1003,1004,1005,1006]),
46+
... y = np.array([2,3,4,4.5,5,5.1]),
47+
... to_exclude = 2,
48+
... )
49+
yhat_u yobs yhat yhat_l polynomial
50+
0 6.817413 5.0 5.500 4.182587 1
51+
1 7.952702 5.1 6.350 4.747298 1
52+
2 9.077182 5.0 4.875 0.672818 2
53+
3 13.252339 5.1 4.975 -3.302339 2
54+
"""
55+
56+
if not isinstance(poly_features, list):
57+
raise ValueError("Please input a list of integers from 0 to 4 for poly_features.")
58+
assert all(0 <= degree <= 4 for degree in poly_features), (
59+
"Please ensure all numbers in poly_features are from 0 to 4."
60+
)
61+
if not isinstance(alpha, float) or 0 > alpha >= 1:
62+
raise ValueError("Please input a float between 0 and 1 for alpha.")
63+
if not isinstance(to_exclude, int):
64+
raise ValueError("Please input an integer between 1 and your sample size for to_exclude.")
65+
assert ((len(t) - to_exclude) >= 4), (
66+
"""The sample size for your model is smaller than 4. This will not produce a good
67+
model. Either reduce to_exclude or increase your sample size to continue."""
68+
)
69+
assert np.isfinite(y).all(), (
70+
f"""Your sample contains missing or infinite values for y at locations
71+
{list(map(tuple, np.where(np.isnan(y))))}. Exclude these values to continue."""
72+
)
73+
assert np.isfinite(t).all(), (
74+
f"""Your sample contains missing or infinite values for t at locations
75+
{list(map(tuple, np.where(np.isnan(t))))}. Exclude these values to continue."""
76+
)
77+
78+
79+
# Sort data by X increasing
80+
idx = np.argsort(t)
81+
t = t[idx]
82+
y = y[idx]
83+
84+
results = pd.DataFrame()
85+
for degree in poly_features:
86+
transforms = make_pipeline(
87+
StandardScaler(),
88+
PolynomialFeatures(degree=degree),
89+
)
90+
91+
# Fit transforms to training data only, apply to all data.
92+
fitted_transforms = transforms.fit(t[:-to_exclude].reshape(-1, 1))
93+
_t = fitted_transforms.transform(t.reshape(-1, 1))
94+
95+
t_train, y_train = _t[:-to_exclude, :], y[:-to_exclude]
96+
t_predict, y_predict = (
97+
_t if predict_all else _t[-to_exclude:, :],
98+
y if predict_all else y[-to_exclude:]
99+
)
100+
101+
# Fit ordinary least squares model to the training data, then predict for the
102+
# prediction data.
103+
model = sm.OLS(y_train, t_train).fit()
104+
yhat = model.predict(t_predict)
105+
106+
# Calculate prediction intervals of fitted model.
107+
_, yhat_l, yhat_u = wls_prediction_std(model, t_predict, alpha=alpha)
108+
109+
# Store model results in master frame
110+
results = results.append(
111+
pd.DataFrame({
112+
"yhat_u" : yhat_u,
113+
"yobs" : y_predict,
114+
"yhat" : yhat,
115+
"yhat_l" : yhat_l,
116+
"polynomial" : degree
117+
}),
118+
ignore_index=True,
119+
)
120+
121+
return results
122+

requirements.txt

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,6 @@
11
numpy>=1.16.0
22
scipy>=0.19.0
33
pandas>=0.24.0
4-
sqlalchemy>=1.3.5
4+
sqlalchemy>=1.3.5
5+
scikit-learn>=0.21.2
6+
statsmodels>=0.10.0

0 commit comments

Comments
 (0)