Рассмотрим ещё один пример построения модели и её оценки с помощью метода наименьших квадратов. Для этого возьмём ещё один встроенный набор данных в R — набор данных по швейцарским кантонам. Поместим в табличку…Давайте табличку на этот раз назовём буковкой t. В табличку t поместим встроенный набор данных swiss. Опять же можно почитать документацию по этому набору данных и обнаружить, что это данные 1888-го года, и переменные здесь: Fertility — это количество детей до пяти лет, делённое на количество женщин, по-моему, до 50-ти лет и отмасштабированное, на 1000 домноженное; Agriculture — это процент мужчин, занятых в сельском хозяйстве; Examination — это какие-то, процент тех, кто получил высокий результат где-то на экзаменах в армии, и так далее. У нас есть описание этого набора данных. И мы попробуем посмотреть, как зависит, соответственно, фертильность, то есть количество детей, деленное на количество женщин, от остальных показателей по кантонам. Опять же, не забывайте никогда смотреть на данные, то есть я рекомендую, во-первых, посмотреть на табличку t вживую. Для этого мы бросим взгляд на неё glimpse(t). Вот мы видим переменные и первые несколько значений. Соответственно, можно посмотреть описательные статистики describe(t). Нажимаем Ctrl Enter, видим среднее значение по каждой переменной. Видим, что у нас 47 наблюдений и так далее. И, конечно, смотрите на графики. К сожалению, когда у вас много количественных переменных, то одной диаграммы рассеяния, которую мы строили, конечно недостаточно. И вообще визуализировать сложный набор данных с большим количеством переменных очень сложно. Но одно из простейших средств посмотреть, как-то хоть себе представить огромное количество переменных — это нарисовать большое количество диаграмм рассеяния. В R есть готовая функция, называется ggpairs, она живёт в пакете GGally, который мы предварительно подгрузили. И можем набрать ggpairs(t) и посмотреть на результат. По графику, например, уже можно понять, что, например, католики не любят протестантов. Вот так вот, раз, и по этому графику вот видите, если я посмотрю переменную доли католического населения в кантоне, то есть резкая поляризация: либо 0 % католического населения, либо 100 %. Ну бывают редкие исключения, но в целом, вот если посмотреть на переменную, то видна поляризация: либо 100 % католического населения, либо 0 %. То есть на этом графике у нас представлена диаграмма рассеяния каждой переменной с каждой: доля католического населения против переменной Education, доля мужчин, занятых в сельском хозяйстве, против переменной «фертильность». И здесь видна некоторая положительная взаимосвязь: чем выше доля мужчин, занятых в сельском хозяйстве, тем выше фертильность, то есть можем сделать вывод, что в сельских регионах, где больше мужчин занято в сельском хозяйстве, там вот выше рождаемость. Соответственно, вот видна какая-то отрицательная взаимосвязь между «фертильностью» и переменной Examination. То есть на этом графике видно много зависимостей, его можно долго изучать. А мы перейдём к оценке модели методом наименьших квадратов. Давайте оценим какую-нибудь модель, пусть будет model2 с помощью метода наименьших квадратов, то есть lm. По данным из таблички t оценим следующую модель: мы посмотрим от чего зависит фертильность. То есть фертильность — это у нас объясняемая переменная. Ставим тильдочку (~) и тут пишем объясняющие переменные. В качестве потенциальных объясняющих возьмём, скажем, долю мужчин, занятых в сельском хозяйстве, переменную education и, скажем, долю католического населения в соответствующем кантоне. Соответственно, запускаем нашу команду. Модель моментально оценивается. Мы можем посмотреть коэффициенты: β с крышкой в оценённой модели. То есть видно, что в нашей модели при прочих равных доля мужчин, занятых в сельском хозяйстве, отрицательно сказывается на фертильности, образование — отрицательно, а доля католического населения положительно сказывается на показателе фертильности. Можно опять же посмотреть спрогнозированные значения фертильности. Можно, соответственно, тут подписаны названия кантонов швейцарских. Можно посмотреть остатки ε с крышкой в нашей model2. Можно посмотреть показатель RSS, deviance(model2). Пропустил буковку n, deviance(model2). И можно посчитать R-квадрат. R-квадрат можно посчитать ещё одним способом. Есть готовое понятие отчёта, которое рассчитывает всякие дополнительные статистики для метода наименьших квадратов. Давайте поместим в список report некий отчёт по модели 2. И теперь из этого отчёта, ну, во-первых, его можно посмотреть. Здесь куча информации, которую мы пока не знаем, но, в частности, среди прочего выводит этот отчёт показатель, который нам уже известен, — множественный R-квадрат 0,64, и его можно из этого списка извлечь, используя $r.squared. Соответственно R-квадрат в этой модели равен 0,64, и он показывает корреляцию, квадрат выборочной корреляции между спрогнозированной переменной и исходной объясняемой переменной. Если я возьму корреляцию из таблички t достану переменную fertility и возьму спрогнозированные значения из model2, это у меня получится корреляция между спрогнозированной переменной и исходной зависимой переменной. И если я возведу эту корреляцию в квадрат, то я получу в точности коэффициент R-квадрат. Аналогично можно в данном наборе данных спрогнозировать: если я возьму новый набор данных nd2, data.frame, куда помещу гипотетические значения для какого-нибудь региона, ну, скажем, рассмотрю регион с долей мужчин, занятых в сельском хозяйстве 0,5, с долей католического населения, тоже пусть будет равной 0,5, и с переменной education, равной 20-ти. Вот создам новый набор данных. И, соответственно, я смогу прогнозировать по model2 для nd2. У меня окажется, что показатель фертильности, то есть это домноженное на 1000 количество детей до пяти лет, делённое на количество женщин. Соответственно, прогноз равен 64.
