import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
%matplotlib inline
from sklearn.model_selection import train_test_split
from sklearn.metrics import mean_squared_error
from sklearn.preprocessing import StandardScaler
from sklearn.model_selection import GridSearchCV
from sklearn.linear_model import LinearRegression
from sklearn.neighbors import KNeighborsRegressor
from sklearn.ensemble import RandomForestRegressor
from sklearn.ensemble import GradientBoostingRegressor
from sklearn.ensemble import AdaBoostRegressor
from sklearn.svm import LinearSVR
from sklearn.learning_curve import learning_curve
from scipy import stats
Данные по успеваемости школьников (средняя школа)
#Нас будет интересовать только успеваемость по математике
#data = pd.read_csv('student-mat.csv')
data = pd.read_csv('../../data/student-mat.csv')
print("Количество учащихся:",len(data))
#Первый взгляд на данные
data.head().T
#В данных пропусков нет, повезло.
print(len(data)-len(data.dropna()))
data.get_ftype_counts()
Всего три типа признаков binary - бинарный, numeric - численный и nominal - категориальный.
school - школа студента (binary: 'GP' - Gabriel Pereira или 'MS' - Mousinho da Silveira)
sex - пол студента (binary: 'F' - female или 'M' - male)
age - возраст студента (numeric: от 15 до 22)
address - тип места проживания (binary: 'U' - город or 'R' - сельская местность)
famsize - размер семьи (binary: 'LE3' - меньше или равно 3 or 'GT3' - больше 3)
Pstatus - живут ли родители совместно (binary: 'T' - живут вместе или 'A' - поотдельности)
Medu - образование матери (numeric: 0 - отсутствует, 1 - начальная школа (4th grade), 2 – с 5-ого по 9ый класс, 3 – оконченная средняя школа or 4 – высшее образование)
Fedu - образование отца (numeric: 0 - отсутствует, 1 - начальная школа (4th grade), 2 – с 5-ого по 9ый класс, 3 – оконченная средняя школа or 4 – высшее образование)
Mjob - работа матери (nominal: 'teacher'-учитель, 'health'-здравоохранения, 'services' - гражданская\государственная работа (например административная или в полиции), 'at_home' - дома или 'other')
Fjob - работа отца (nominal: 'teacher'-учитель, 'health'-здравоохранения, 'services' - гражданская\государственная работа (например административная или в полиции), 'at_home' - дома или 'other')
reason - причина выбора школы (nominal:'home'-близко к дому, 'reputation'-репутация школы, 'course'-предпочтения в курсах или 'other')
guardian - представитель студента (nominal: 'mother', 'father' или 'other')
traveltime - время от дома до школы (numeric: 1 - меньше 15 мин., 2 - от 15 до 30 мин., 3 - от 30 мин. до 1 часа, или 4 - больше 1 часа)
studytime - еженедельные временные затраты на обучение (numeric: 1 - меньше 2 часов, 2 - от 2 до 5 часов, 3 - от 5 to 10 hours, или 4 - больше 10 часов)
failures - количество проваленных классов в прошлом (numeric: n если n=1,2,3, иначе 4)
schoolsup - дополнительная поддержа в обучении (видимо финансовая) (binary: yes or no)
famsup - семейная образовательная поддержка (binary: yes or no)
paid - оплата дополнительных уроков по предмету (Math) (binary: yes or no)
activities - внеклассовые активности (binary: yes or no)
nursery - учился ли в детском саду (binary: yes or no)
higher - желание получить высшее образование (binary: yes or no)
internet - домашний доступ к интернету (binary: yes or no)
romantic - в любовных отношениях (binary: yes or no)
famrel - качество семейных взаимоотношений (numeric: от 1 - очень плохие до 5 - отличные)
freetime - кол-во свободного времени после школы (numeric: от 1 - очень мало до 5 - очень много)
goout - как часто гуляет с друзьями (numeric: от 1 - очень мало до 5 - очень много)
Dalc - кол-во потребляемого алкоголя в будни (numeric: от 1 - очень мало до 5 - очень много)
Walc - кол-во потребляемого алкоголя в выходные (numeric: от 1 - очень мало до 5 - очень много)
health - состояние здоровья (numeric: от 1 - очень плохое до 5 - отличное)
absences - количество пропусков занятий (numeric: от 0 до 93)
Семестровые оценки по математике (целевые признаки):
G1 - оценка за первый семестр (numeric: от 0 до 20)
G2 - оценка за второй семестр (numeric: от 0 до 20)
G3 - финальная оценка (numeric: от 0 до 20)
Всего признаков:33
Количество учащихся:395
Пожалуй, интересней всего в этой задаче посмотреть взаимодействия между признаками и выделить признаки наиболее влияющие на успеваемость. Научится предсказывать оценку тоже не плохо, например, зная оценки в одном классе школы можно предсказать оценки в другом классе.
binary_features=['school','sex','address','famsize','Pstatus','schoolsup','famsup','paid','activities','nursery','higher','internet','romantic']
nominal_features=['Mjob','Fjob','reason','guardian']
numeric_features=['age','Medu','Fedu','traveltime','studytime','failures','famrel','freetime','goout','Dalc','Walc','health','absences']
target_features=['G1','G2','G3']
len(binary_features)+len(nominal_features)+len(numeric_features)+len(target_features)
data_num=data[numeric_features]
data_bin=data[binary_features]
data_nom=data[nominal_features]
data_target=data[target_features]
data_bin_dum=pd.get_dummies(data_bin)
#Бинарные признаки разбились на взаимодополняющие пары, из каждой пары достаточно оставить только один признак.
data_bin_cut=data_bin_dum.iloc[:,1::2]
data_nom_dum=pd.get_dummies(data_nom)
#Зная другие признаки соответствующие исходному категориальному признаку, значения этих признаков определяется однозначно.
#Удаляем
data_nom_cut=data_nom_dum.drop(['Mjob_other','Fjob_other','reason_other','guardian_other'],axis=1)
data_reformed=pd.concat([data_num,data_bin_cut,data_nom_cut,data_target],axis=1)
#Целевая переменная, которую мы будем предсказывать - это средняя сумма баллов за все экзамены.
data_reformed['G_average']=(data_reformed['G1']+data_reformed['G2']+data_reformed['G3'])/3
#Посмотрим на статистику
data_reformed.describe().T
#В бинарных признаках например: higher_yes 1 - означает да(хочет получить высшее образование), 0 - означает нет.
sns.set(font_scale=2.5)
plt.figure(figsize=(40,30))
corr_matrix=data_reformed.corr()
sns.heatmap(corr_matrix,annot=True,fmt = ".2f",cbar = True,cmap='PuOr',annot_kws={"size":18})
#Присутствует также сильная корреляция между уровнем образования матери и отца
print(corr_matrix['Fedu']['Medu'])
#Есть и вполне естественные корреляции между потреблением алкоголя по будням, потреблением алкоголя в выходные и тем
#сколько человек гуляет с друзьями.
print(corr_matrix['Dalc']['Walc'])
print(corr_matrix['Walc']['goout'])
# Очень интересно, большая отрицательная корреляция между мужским полом и временем затрат на учебу
print(corr_matrix['sex_M']['studytime'])
# А также вполне логичная корреляция между проживанием в городе (а не сельской местности) и временем на дорогу до школы
print(corr_matrix['address_U']['traveltime'])
#Бросается в глаза сильная корелляция между количеством не сдач экзамена ранее и средней оценкой
print(corr_matrix['G_average']['failures'])
#В некотором смысле failures это обратный признак к средней оценке, пожалуй, для чистоты эксперимента его надо исключить
data_reformed=data_reformed.drop(['failures'],axis=1)
#Сильная положительная корреляция оценки с уровнем образования матери и отца, временем занятий и желанием получить
#высшее образование.
print(corr_matrix['G_average']['Medu'])
print(corr_matrix['G_average']['Fedu'])
print(corr_matrix['G_average']['studytime'])
print(corr_matrix['G_average']['higher_yes'])
#А также оценка по математике положительно коррелирует с мужским полом.
print(corr_matrix['sex_M']['G_average'])
#Из отрицательных корреляций наиболее выделяются тусовки с друзьями, schoolsup и работа матери на дому, ого!
print(corr_matrix['G_average']['goout'])
print(corr_matrix['G_average']['schoolsup_yes'])
print(corr_matrix['G_average']['Mjob_at_home'])
#Есть подозрение, что schoolsup - это поддержка бедных детей (а бедные дети скорее всего хуже учатся).
data_reformed.columns
features_to_hist=['age','Medu','Fedu','studytime']
#Строим графики
sns.set(font_scale=1.2)
data_reformed[features_to_hist].hist(figsize=(12,10),bins=8)
#sns.pairplot(data_reformed[features_to_hist + ['G_average']])
#Много родителей не закончили среднюю школу, гм!
#Большинство учеников тратят на учебу от 2 до 5 часов в неделю
#А средний возраст равен:
print('средний возраст:',data_reformed['age'].mean())
#Людей старше 19 всего лишь 5, максимальный возраст - 22 года.
data_reformed['age'].value_counts()
#Существенных выбросов в данных тоже не обнаружено
data_reformed['G_average'].hist(bins=25)
# Немного напоминает нормальное распределение
# Посмотрим, что говорят статистические тесты.
print(stats.normaltest(data_reformed['G_average']))
print('skew=',stats.skew(data_reformed['G_average']))
print(stats.skewtest(data_reformed['G_average']))
#Отделение признаков от целевой переменной
X=data_reformed.iloc[:,:-4]
y=data_reformed.iloc[:,-1:]
#Заодно сразу сделаем разбиение на train и test.
#Обязательно перемешать!
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.30,shuffle=True,random_state=42)
#Среднее значение оценки на трейне
y_mean=y_train.mean()
y_mean_for_test=[float(y_mean) for x in range(len(y_test))]
y_mean_for_train=[float(y_mean) for x in range(len(y_train))]
#Результаты приближения средним на трейне
print(mean_squared_error(y_train,y_mean_for_train))
#Результаты приближения средним на тесте
baseline_mean=mean_squared_error(y_test,y_mean_for_test)
print(mean_squared_error(y_test,y_mean_for_test))
#Для удобства сделаем массивы целевых переменных строками
y_train_row=[]
y_test_row=[]
for x in y_train.iloc[:,0]:
y_train_row.append(x)
for x in y_test.iloc[:,0]:
y_test_row.append(x)
#Обучаем базовый случайны лес для отбора признаков
rfc_base=RandomForestRegressor(n_estimators=500,random_state=42)
rfc_base.fit(X_train,y_train_row)
#Выводим значимость признаков
features = pd.DataFrame(rfc_base.feature_importances_, index=X_train.columns,
columns=['Importance']).sort_values(['Importance'], ascending=False)
features
#Отбросим 15 самых малозначимых признаков
#Я пробовал удалять разное количество малозначимых признаков и удалить 15 оказалось оптимальным.(с точки зрения кроссвалидации, см. далее)
features_cutted=features.iloc[:-15].index
features_cutted
#Оставим в трейне и тесте только эти признаки
X_train_cutted=X_train[features_cutted]
X_test_cutted=X_test[features_cutted]
#На этом преобработка данных полностью закончена
#Приступим к выбору модели
#Список регрессоров
regressors = [LinearRegression(),
GradientBoostingRegressor(random_state=42),
RandomForestRegressor(random_state=42),
LinearSVR(random_state=42)]
regressor_name = ['LinearRegression',
'GradientBoostingRegressor',
'RandomForestRegressor',
'LinearSVR']
#Параметры к регрессорам
scores = []
fits = []
linear_params = {'normalize': (True, False)}
gbr_params = {'n_estimators': [100, 300, 500],
'learning_rate':(0.1, 0.5, 1),
'max_depth': list(range(3, 10, 2)),
'min_samples_leaf': list(range(10, 31, 10))}
forest_params = {'n_estimators': [100, 300, 500],
'max_depth': list(range(3, 10, 2)),
'min_samples_leaf': list(range(10, 31, 10))}
svm_params = {'loss' : ('epsilon_insensitive', 'squared_epsilon_insensitive'), 'C': (.5, 1, 2)}
params = [linear_params, gbr_params, forest_params, svm_params]
#Перебираем параметры регрессоров в поисках лучшего (на 5 фолдах)
np.random.seed(42)
for i, each_regressor in enumerate(regressors):
reg = each_regressor
reg_params = params[i]
grid = GridSearchCV(reg, reg_params,
cv=5,
scoring='neg_mean_squared_error',
n_jobs=-1)
grid.fit(X_train_cutted, y_train_row)
fits.append(grid.best_params_)
reg_best_score = grid.best_score_
scores.append(reg_best_score)
print(regressor_name[i], -reg_best_score, "\n", grid.best_params_, "\n")
LinearRegression 14.0748985397 {'normalize': False}
GradientBoostingRegressor 11.7584062474 {'learning_rate': 0.1, 'max_depth': 3, 'min_samples_leaf': 20, 'n_estimators': 100}
RandomForestRegressor 11.3097947653 {'max_depth': 9, 'min_samples_leaf': 10, 'n_estimators': 100}
LinearSVR 12.8721586628 {'C': 2, 'loss': 'epsilon_insensitive'}
#Углубленный подбор гиперпараметров для случайного леса
np.random.seed(42)
forest_params_deep = {'n_estimators': [100,150,200,300,500], #n_estimators - количество деревьев в случайном лесе
'max_depth': list(range(3, 13, 2)), #max_depth - максимальная глубина дерева
'min_samples_leaf': list(range(5, 30, 5))}#min_samples_leaf - минимальное количество объектов в листе дерева.
rfr=RandomForestRegressor(random_state=42)
grid_rfr = GridSearchCV(rfr, forest_params_deep,
cv=5,
scoring='neg_mean_squared_error',
n_jobs=-1)
grid_rfr.fit(X_train_cutted, y_train_row)
#Результат стал еще лучше!
print(-grid_rfr.best_score_,'\n',grid_rfr.best_params_)
#Еще раз взглянем на важность признаков и держа в голове ранее исследованную матрицу корреляции, построим новые признаки.
features = pd.DataFrame(grid_rfr.best_estimator_.feature_importances_, index=X_train_cutted.columns,
columns=['Importance']).sort_values(['Importance'], ascending=False)
features
X_train_new_features=X_train_cutted.copy()
X_test_new_features=X_test_cutted.copy()
#Вспомним про сильно коррелированные признаки потребления алкоголя в будни и выходные, пожалуй их стоит объединить в один.
X_train_new_features['Talc']=X_train_new_features['Dalc']+X_train_new_features['Walc']
X_test_new_features['Talc']=X_test_new_features['Dalc']+X_test_new_features['Walc']
X_train_new_features=X_train_new_features.drop(['Dalc','Walc'],axis=1)
X_test_new_features=X_test_new_features.drop(['Dalc','Walc'],axis=1)
#Образование матери и отца тоже сильно коррелирует, пожалуй стоит учитывать только суммарное образование родителей.
X_train_new_features['Pedu']=X_train_new_features['Medu']+X_train_new_features['Fedu']
X_test_new_features['Pedu']=X_test_new_features['Medu']+X_test_new_features['Fedu']
X_train_new_features=X_train_new_features.drop(['Medu','Fedu'],axis=1)
X_test_new_features=X_test_new_features.drop(['Medu','Fedu'],axis=1)
#Вспомним, что потребление алкоголя кореллирует с временем прогулок с друзьями, попробуем перемножить эти признаки.
X_train_new_features['goout_alc']=X_train_new_features['goout']*X_train_new_features['Talc']
X_test_new_features['goout_alc']=X_test_new_features['goout']*X_test_new_features['Talc']
X_train_new_features=X_train_new_features.drop(['goout','Talc'],axis=1)
X_test_new_features=X_test_new_features.drop(['goout','Talc'],axis=1)
#Есть подозрение, что так как мужской пол негативно кореллирует с затратами времени на учебу, но мужской пол
#положительно кореллирует с оценками по математике, то возможно
#что если мальчик тратит больше времени на учебу, то его результат усиливается сильнее чем у девочек.
X_train_new_features['studytime_eff']=X_train_new_features['studytime']*(X_train_new_features['sex_M']*0.5 + 1)
X_test_new_features['studytime_eff']=X_test_new_features['studytime']*(X_test_new_features['sex_M']*0.5 + 1)
X_train_new_features=X_train_new_features.drop(['studytime','sex_M'],axis=1)
X_test_new_features=X_test_new_features.drop(['studytime','sex_M'],axis=1)
X_train_new_features.T
#Стало ли лучше? Снова проведем gridsearchCV.
np.random.seed(42)
forest_params_deep = {'n_estimators': [100,150,200,300,500],
'max_depth': list(range(3, 13, 2)),
'min_samples_leaf': list(range(5, 30, 5))}
rfr=RandomForestRegressor(random_state=42)
grid_rfr = GridSearchCV(rfr, forest_params_deep,
cv=5,
scoring='neg_mean_squared_error',
n_jobs=-1)
grid_rfr.fit(X_train_new_features, y_train_row)
#И результат снова немного улучшается
print(-grid_rfr.best_score_,'\n',grid_rfr.best_params_)
#Снова смотрим на важность признаков
features = pd.DataFrame(grid_rfr.best_estimator_.feature_importances_, index=X_train_new_features.columns,
columns=['Importance']).sort_values(['Importance'], ascending=False)
features
def plot_with_std(x, data, **kwargs):
mu, std = data.mean(1), data.std(1)
lines = plt.plot(x, mu, '-', **kwargs)
plt.fill_between(x, mu - std, mu + std, edgecolor='none',
facecolor=lines[0].get_color(), alpha=0.2)
def plot_learning_curve(reg, X, y, scoring, cv=5):
train_sizes = np.linspace(0.05, 1, 20)
n_train, val_train, val_test = learning_curve(reg,
X, y, train_sizes, cv=cv,
scoring=scoring)
plot_with_std(n_train, val_train, label='training scores', c='green')
plot_with_std(n_train, val_test, label='validation scores', c='red')
plt.xlabel('Training Set Size'); plt.ylabel(scoring)
plt.legend()
#Это график MSE со знаком минус, поэтому функции возрастают.
plot_learning_curve(grid_rfr.best_estimator_,
X_train_new_features, y_train_row, scoring='neg_mean_squared_error', cv=5)
#Впринципе результаты сравнимы с результатами на кроссвалидации.
#Хотя все же немного обидно, что различие существенно, ведь разбиение на трейн и тест было случайным.
#Исходное приближение средним улучшить удалось.
print(mean_squared_error(y_test_row,grid_rfr.best_estimator_.predict(X_test_new_features)),'vs baseline:',baseline_mean)
Нам удалось выделить наиболее значимые признаки влияющие на успеваемость школьников и это, наверное, самое главное, потому что на признаки можно повлиять. Также, зная успеваемость и признаки школьников из одного класса, можно спрогнозировать успеваемость в другом классе. Интересно было бы исследовать важность признаков для разных школ. Например, какие признаки наиболее влияют на успеваемость для сельской школы, а какие для городской. Для улучшения решения можно было бы дополнительно использовать градиентный бустинг, который себя неплохо показал, а потом его результаты объединить со случайным лесом.