66from 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