Тест. Критерии Стьюдента


In [1]:
from __future__ import division

import numpy as np
import pandas as pd

from sklearn.model_selection import train_test_split
from sklearn.linear_model import LinearRegression
from sklearn.ensemble import RandomForestRegressor

from scipy import stats
from statsmodels.stats.weightstats import CompareMeans, DescrStatsW

from IPython.core.interactiveshell import InteractiveShell
InteractiveShell.ast_node_interactivity = "all"

Уровень кальция в крови здоровых молодых женщин равен в среднем 9.5 милиграммам на децилитр и имеет характерное стандартное отклонение 0.4 мг/дл. В сельской больнице Гватемалы для 160 здоровых беременных женщин при первом обращении для ведения беременности был измерен уровень кальция; среднее значение составило 9.57 мг/дл. Можно ли утверждать, что средний уровень кальция в этой популяции отличается от 9.5?

Посчитайте достигаемый уровень значимости. Поскольку известны только среднее и дисперсия, а не сама выборка, нельзя использовать стандартные функции критериев — нужно реализовать формулу достигаемого уровня значимости самостоятельно.

Округлите ответ до четырёх знаков после десятичной точки.


In [2]:
exp_calc = 9.5
st_dev_calc = 0.4

In [3]:
num = 160
mean_calc = 9.57

In [4]:
def z_test(mean_val, exp_val, st_dev, num):
    standard_error = st_dev / np.sqrt(num)
    return (mean_val - exp_val) / standard_error

In [5]:
z_stat_calc = z_test(mean_calc, exp_calc, st_dev_calc, num)
z_stat_calc


Out[5]:
2.2135943621178749

In [6]:
def p_val_double(z_stat):
    return 2 * (1 - stats.norm.cdf(z_stat))

In [7]:
print('p-value: %.4f' % p_val_double(z_stat_calc))


p-value: 0.0269

Имеются данные о стоимости и размерах 53940 бриллиантов:

diamonds.txt

Отделите 25% случайных наблюдений в тестовую выборку с помощью функции sklearn.cross_validation.train_test_split (зафиксируйте random state = 1). На обучающей выборке настройте две регрессионные модели:

линейную регрессию с помощью LinearRegression без параметров случайный лес с помощью RandomForestRegressor с random_state=1. Какая из моделей лучше предсказывает цену бриллиантов? Сделайте предсказания на тестовой выборке, посчитайте модули отклонений предсказаний от истинных цен. Проверьте гипотезу об одинаковом среднем качестве предсказаний, вычислите достигаемый уровень значимости. Отвергается ли гипотеза об одинаковом качестве моделей против двусторонней альтернативы на уровне значимости α=0.05?


In [8]:
diamonds = pd.read_table('diamonds.txt')
diamonds.describe()
diamonds.head()
diamonds.info()


Out[8]:
carat depth table price x y z
count 53940.000000 53940.000000 53940.000000 53940.000000 53940.000000 53940.000000 53940.000000
mean 0.797940 61.749405 57.457184 3932.799722 5.731157 5.734526 3.538734
std 0.474011 1.432621 2.234491 3989.439738 1.121761 1.142135 0.705699
min 0.200000 43.000000 43.000000 326.000000 0.000000 0.000000 0.000000
25% 0.400000 61.000000 56.000000 950.000000 4.710000 4.720000 2.910000
50% 0.700000 61.800000 57.000000 2401.000000 5.700000 5.710000 3.530000
75% 1.040000 62.500000 59.000000 5324.250000 6.540000 6.540000 4.040000
max 5.010000 79.000000 95.000000 18823.000000 10.740000 58.900000 31.800000
Out[8]:
carat depth table price x y z
0 0.23 61.5 55.0 326 3.95 3.98 2.43
1 0.21 59.8 61.0 326 3.89 3.84 2.31
2 0.23 56.9 65.0 327 4.05 4.07 2.31
3 0.29 62.4 58.0 334 4.20 4.23 2.63
4 0.31 63.3 58.0 335 4.34 4.35 2.75
<class 'pandas.core.frame.DataFrame'>
RangeIndex: 53940 entries, 0 to 53939
Data columns (total 7 columns):
carat    53940 non-null float64
depth    53940 non-null float64
table    53940 non-null float64
price    53940 non-null int64
x        53940 non-null float64
y        53940 non-null float64
z        53940 non-null float64
dtypes: float64(6), int64(1)
memory usage: 2.9 MB

In [9]:
diamonds.columns


Out[9]:
Index([u'carat', u'depth', u'table', u'price', u'x', u'y', u'z'], dtype='object')

In [10]:
X_diam = diamonds.drop(['price'], axis=1)
X_diam.shape
y_diam = diamonds.loc[:, diamonds.columns == 'price']
np.ravel(y_diam).shape


Out[10]:
(53940, 6)
Out[10]:
(53940L,)

In [11]:
X_diam_train, X_diam_test, y_diam_train, y_diam_test = train_test_split(X_diam, y_diam, random_state=1)

In [12]:
clf_lr = LinearRegression()
clf_lr.fit(X_diam_train, y_diam_train)


Out[12]:
LinearRegression(copy_X=True, fit_intercept=True, n_jobs=1, normalize=False)

In [13]:
pred_price_lr = clf_lr.predict(X_diam_test)

In [14]:
pred_price_lr_sub = pred_price_lr - y_diam_test
pred_price_lr_sub.describe()


Out[14]:
price
count 13485.000000
mean 19.229235
std 1463.058136
min -12455.940789
25% -342.670547
50% 63.649682
75% 652.518106
max 18239.846360

In [15]:
clf_rf = RandomForestRegressor(random_state=1)
clf_rf.fit(X_diam_train, y_diam_train.values.ravel())


Out[15]:
RandomForestRegressor(bootstrap=True, criterion='mse', max_depth=None,
           max_features='auto', max_leaf_nodes=None,
           min_impurity_split=1e-07, min_samples_leaf=1,
           min_samples_split=2, min_weight_fraction_leaf=0.0,
           n_estimators=10, n_jobs=1, oob_score=False, random_state=1,
           verbose=0, warm_start=False)

In [16]:
pred_price_rf = clf_rf.predict(X_diam_test)

In [17]:
pred_price_rf_sub = pred_price_rf[:, np.newaxis] - y_diam_test
pred_price_rf_sub.describe()


Out[17]:
price
count 13485.000000
mean 48.200849
std 1407.659287
min -12710.200000
25% -245.200000
50% 33.700000
75% 402.100000
max 9124.200000

In [21]:
stats.ttest_rel(np.abs(pred_price_lr_sub), np.abs(pred_price_rf_sub))


Out[21]:
Ttest_relResult(statistic=array([ 12.74505678]), pvalue=array([  5.42865482e-37]))

In [25]:
print "95%% confidence interval: [%f, %f]" % DescrStatsW(np.abs(pred_price_lr_sub) - np.abs(pred_price_rf_sub)).tconfint_mean()


95% confidence interval: [72.497266, 98.849770]