做了一段時(shí)間分析后兄世，你是否會(huì)感覺分析結(jié)果還停留在數(shù)據(jù)表面現(xiàn)象，如果你和我有同感邻薯，請耐心讀完本篇文章裙戏，會(huì)受益匪淺。接下來一段時(shí)間，本人也會(huì)深入研究統(tǒng)計(jì)學(xué)方面的知識(shí)。因?yàn)閷W(xué)了統(tǒng)計(jì)學(xué)二跋，你會(huì)發(fā)現(xiàn)很多時(shí)候的分析并不靠譜。比如很多人都喜歡用平均數(shù)去分析一個(gè)事物的結(jié)果信柿，但是這往往是粗糙的，不準(zhǔn)確的醒第。如果學(xué)了統(tǒng)計(jì)學(xué)渔嚷，那么我們就能以更多更科學(xué)的角度看待數(shù)據(jù)。更深層次的挖掘出數(shù)據(jù)背后的價(jià)值稠曼！統(tǒng)計(jì)學(xué)是數(shù)據(jù)分析的基石形病。
本文通過使用Python+Pandas+Statsmodels建立簡單一元線性回歸模型、多元線性回歸模型來預(yù)測世界幸福指數(shù)霞幅。
主要內(nèi)容源自博客：https://www.sohu.com/a/419787355_568359
我在原文基礎(chǔ)上增加了大量知識(shí)點(diǎn)的說明及解釋漠吻，以達(dá)到學(xué)習(xí)線性回歸的效果。最后的預(yù)測結(jié)果因?yàn)閿?shù)據(jù)集的局限司恳，結(jié)果僅供參考途乃，主要目的還是通過一個(gè)實(shí)際的案例，了解使用線性回歸建模的過程扔傅，測試耍共、分析。
一猎塞、項(xiàng)目背景
《世界幸福指數(shù)報(bào)告》是對(duì)全球幸福狀況的一次具有里程碑意義的調(diào)查试读。
民意測驗(yàn)機(jī)構(gòu)蓋洛普從2012年起，每年都會(huì)在聯(lián)合國計(jì)劃下發(fā)布《世界幸福指數(shù)報(bào)告》荠耽，報(bào)告會(huì)綜合兩年內(nèi)150多個(gè)國家的國民對(duì)其所處社會(huì)钩骇、城市和自然環(huán)境等因素進(jìn)行評(píng)價(jià)后，再根據(jù)他們所感知的幸福程度對(duì)國家進(jìn)行排名铝量。
《世界幸福指數(shù)報(bào)告》的編撰主要依賴于對(duì)150多個(gè)國家的1000多人提出一個(gè)簡單的主觀性問題：“如果有一個(gè)從0分到10分的階梯伊履，頂層的10分代表你可能得到的最佳生活，底層的0分代表你可能得到的最差生活款违。你覺得你現(xiàn)在在哪一層？”
那么哪個(gè)國家在總體幸福指數(shù)上排名最高群凶？哪些因素對(duì)幸福指數(shù)的影響最大插爹？今天我們就用Python+Pandas+Statsmodels來聊一聊。
二、數(shù)據(jù)集說明
- rank：幸福指數(shù)排名
- region：國家
- happiness：幸福指數(shù)得分
- gdp_per_capita：GDP（人均國內(nèi)生產(chǎn)總值）
- healthy_life_expectancy：健康預(yù)期壽命
- freedom_to_life_choise：自由權(quán)
- generosity：慷慨程度
- year：年份
- corruption_perceptions：清廉指數(shù)
- social_support：社會(huì)支持（客觀上物質(zhì)上的援助和直接服務(wù)赠尾；主觀上指個(gè)體感到在社會(huì)中被尊重力穗、被支持和被理解的情緒體驗(yàn)和滿意程度。）
二气嫁、什么是線性回歸当窗？
線性回歸是利用數(shù)理統(tǒng)計(jì)中回歸分析，來確定兩種或兩種以上變量間相互依賴的定量關(guān)系的一種統(tǒng)計(jì)分析方法寸宵。
三崖面、變量選擇
【因變量】 我們選擇happiness(幸福指數(shù))作為因變量，該變量度量了各個(gè)國家的幸福指數(shù)梯影。
【自變量】我們選擇了影響幸福的六個(gè)因素巫员，比如GDP、健康預(yù)期壽命甲棍、自由權(quán)简识、慷慨程度、清廉指數(shù)感猛、社會(huì)支持七扰。
四、 數(shù)據(jù)導(dǎo)入和數(shù)據(jù)整理
import pandas as pd
import numpy as np
# 可視化
import matplotlib.pyplot as plt
import seaborn as sns
import plotly as py
import plotly.graph_objs as go
import plotly.express as px
from plotly.offline import init_notebook_mode, iplot, plot
init_notebook_mode(connected=True)
plt.style.use('seaborn')
# 讀入數(shù)據(jù)
df_2015 = pd.read_csv('./deal_data/2015.csv')
df_2016 = pd.read_csv('./deal_data/2016.csv')
df_2017 = pd.read_csv('./deal_data/2017.csv')
df_2018 = pd.read_csv('./deal_data/2018.csv')
df_2019 = pd.read_csv('./deal_data/2019.csv')
# 新增列-年份
df_2015["year"] = str(2015)
df_2016["year"] = str(2016)
df_2017["year"] = str(2017)
df_2018["year"] = str(2018)
df_2019["year"] = str(2019)
# 合并數(shù)據(jù)
df_all = df_2015.append([df_2016, df_2017, df_2018, df_2019], sort=False)
df_all.drop('Unnamed: 0', axis=1, inplace=True)
df_all.head()
print(df_2015.shape, df_2016.shape, df_2017.shape, df_2018.shape, df_2019.shape)
五陪白、數(shù)據(jù)可視化
1颈走、2019世界幸福地圖
data = dict(type = 'choropleth',
locations = df_2019['region'],
locationmode = 'country names',
colorscale = 'RdYlGn',
z = df_2019['happiness'],
text = df_2019['region'],
colorbar = {'title':'Happiness'})
layout = dict(title = 'Geographical Visualization of Happiness Score in 2019',
geo = dict(showframe = True, projection = {'type': 'azimuthal equal area'}))
choromap3 = go.Figure(data = [data], layout=layout)
plot(choromap3, filename='./html/世界幸福地圖.html')
結(jié)論：整體來看，北歐的國家幸福指數(shù)較高拷泽，如冰島疫鹊、丹麥、挪威司致、芬蘭拆吆；東非和西非的國家幸福指數(shù)較低，如多哥脂矫、布隆迪枣耀、盧旺達(dá)和坦桑尼亞。
2庭再、2019世界幸福國家排行Top10
# 合并數(shù)據(jù)
rank_top10 = df_2019.head(10)[['rank', 'region', 'happiness']]
last_top10 = df_2019.tail(10)[['rank', 'region', 'happiness']]
rank_concat = pd.concat([rank_top10, last_top10])
# 條形圖
fig = px.bar(rank_concat,
x="region",
y="happiness",
color="region",
title="World's happiest and least happy countries in 2019")
plot(fig, filename='./html/2019世界幸福國家排行Top10和Last10.html')
3捞奕、幸福指數(shù)相關(guān)性
# 熱力圖
plt.figure(figsize=(25, 20))
sns.heatmap(df_all.corr(), cmap='rainbow', linewidths=0.1, annot=True)
plt.title('Correlation between numeric variables', fontsize=18)
plt.xticks(fontsize=13)
plt.yticks(fontsize=13)
plt.show()
結(jié)論：從影響因素相關(guān)性熱力圖可以看出，在影響幸福得分的因素中拄轻，GDP颅围、社會(huì)支持、健康預(yù)期壽命呈現(xiàn)高度相關(guān)恨搓，自由權(quán)呈現(xiàn)中度相關(guān)院促，國家的廉政水平呈現(xiàn)低度相關(guān)筏养，慷慨程度則呈現(xiàn)極低的相關(guān)性；
GDP與健康預(yù)期壽命常拓、社會(huì)支持之間存在高度相關(guān)渐溶。說明GDP高的國家，醫(yī)療水平和社會(huì)福利較為完善弄抬，人民的預(yù)期壽命也會(huì)越高茎辐；
健康預(yù)期壽命與社會(huì)支持之間存在中度相關(guān)性。
4掂恕、以下分別觀察各個(gè)因素的影響程度
①GDP和幸福得分
# 散點(diǎn)圖
fig = px.scatter(df_all, x='gdp_per_capita',
y='happiness',
facet_row='year',
color='year',
trendline='ols'
)
fig.update_layout(height=800, title_text='GDP per capita and Happiness Score')
plot(fig, filename='./html/GDP和幸福得分.html')
結(jié)論：人均GDP與幸福得分呈高度線性正相關(guān)關(guān)系拖陆，GDP越高的國家，幸福水平相對(duì)越高
②健康預(yù)期壽命和幸福得分
# 散點(diǎn)圖
fig = px.scatter(df_all, x='healthy_life_expectancy',
y='happiness',
facet_row='year',
color='year',
trendline='ols'
)
fig.update_layout(height=800, title_text='Healthy Life Expecancy and Happiness Score')
plot(fig, filename='./html/健康預(yù)期壽命和幸福得分.html')
結(jié)論：健康預(yù)期壽命與幸福得分呈高度線性正相關(guān)關(guān)系竹海，健康預(yù)期壽命越高的國家慕蔚，幸福水平相對(duì)越高。
③動(dòng)態(tài)圖展示（GDP&happiness）
fig = px.scatter(df_all,
x='healthy_life_expectancy',
y='happiness',
animation_frame='year',
animation_group='region',
size='rank',
color='region',
hover_name='region',
trendline='ols'
)
fig.update_layout(title_text='Happiness Rank vs healthy_life_expectancy')
plot(fig, filename='./html/健康預(yù)期壽命和幸福水平動(dòng)態(tài)圖展示.html')
④動(dòng)態(tài)圖展示（healthy_life_expectancy&happiness）
fig = px.scatter(df_all,
x='healthy_life_expectancy',
y='happiness',
animation_frame='year',
animation_group='region',
size='rank',
color='region',
hover_name='region',
trendline='ols'
)
fig.update_layout(title_text='Happiness Rank vs healthy_life_expectancy')
plot(fig, filename='./html/健康預(yù)期壽命和幸福水平動(dòng)態(tài)圖展示.html')
五斋配、數(shù)據(jù)建模
我們使用線性回歸進(jìn)行建立一個(gè)基準(zhǔn)模型孔飒，首先篩選一下建模變量，并刪除空值記錄艰争。
from statsmodels.formula.api import ols
sel_cols = ['happiness', 'gdp_per_capita', 'healthy_life_expectancy',
'freedom_to_life_choise', 'corruption_perceptions', 'generosity']
# 重置索引
df_model.index = range(df_model.shape[0])
df_model = df_all[sel_cols]
# 刪除空值
df_model = df_model.dropna()
df_model.head()
1坏瞄、ols最小二乘法
本文使用最小二乘法來建模，它是很基礎(chǔ)甩卓、很強(qiáng)大的方法鸠匀，應(yīng)用非常廣泛。
最小二乘法通過使得因變量觀測值與因變量估計(jì)值之間的殘差平方和達(dá)到最小的方法逾柿，得到估計(jì)的回歸方程缀棍。最小二乘準(zhǔn)則就是選擇能與樣本數(shù)據(jù)有最佳擬合方程的準(zhǔn)則。
最小二乘法建立在假定之上机错，如果假設(shè)成立爬范，建立的模型可以用來較準(zhǔn)確地預(yù)測數(shù)據(jù)。反之弱匪，假設(shè)不成立青瀑，模型論斷將失去其有效性。
最小二乘法假定：
1）線性關(guān)系：假定因變量與自變量之間存在線性關(guān)系萧诫。如果不存在線性關(guān)系斥难，線性回歸不是解釋數(shù)據(jù)的正確模型。
2）無多重共線性：因變量之間不存在相關(guān)性帘饶。如果某些因變量之間存在密切聯(lián)系哑诊，可以嘗試刪去其中一個(gè)或多個(gè)相關(guān)的因變量。因?yàn)槎嘤嗟囊蜃兞刻峁┝巳哂嘈畔⒓翱蹋蕹嘤嘧兞坎⒉粫?huì)大大降低修正判定系數(shù)搭儒。
3）零條件均值假定：觀測值和線性擬合估計(jì)值之間的平均殘差為0穷当。有時(shí)候觀測值和估計(jì)值相比，偏大淹禾；有時(shí)候偏小，但之間的殘差不會(huì)無控制地偏向于一群值茴扁。
4）同方差：對(duì)于自變量的不同取值铃岔，因變量的誤差項(xiàng)都是獨(dú)立的，方差是相同的峭火。
5）無自相關(guān)性/序列相關(guān)：自相關(guān)是指一個(gè)變量同自身其他觀測值有相互關(guān)系毁习。比如，如果今天的股票價(jià)格影響著明天的股票價(jià)格卖丸，那么股票價(jià)格就是序列相關(guān)的纺且。
2、簡單線性回歸
回歸模型
簡單線性回歸使用一個(gè)自變量來預(yù)測一個(gè)因變量稍浆，二者之間的關(guān)系可以用一條直線近似表示载碌。
簡單線性回歸模型：
建模
我們將使用statsmodels中ols功能，構(gòu)建happiness同gdp_per_capita之間的模型衅枫。
Statsmodels是一個(gè)很強(qiáng)大的Python庫嫁艇，用于擬合多種統(tǒng)計(jì)模型，執(zhí)行統(tǒng)計(jì)測試以及數(shù)據(jù)探索和可視化弦撩。
對(duì)于線性回歸linear regression步咪，我們可以使用Statsmodels庫中最小二乘法OLS（Ordinary-Least-Square）的功能來實(shí)現(xiàn)，可以得到豐富的數(shù)據(jù)信息益楼。
# 建立簡單線性回歸模型
lm_m = ols(formula='happiness ~ gdp_per_capita',
data=df_model).fit()
lm_m.summary()
1） 修正判定系統(tǒng)Adj.R-squared：63.3%猾漫。幸福指數(shù)變異性的63.3%能被其與gdp_per_capita之間的線性關(guān)系解釋。
2）回歸系數(shù)：2.19感凤。代表GDP每增加一個(gè)單位悯周，幸福指數(shù)將增加2.19個(gè)單位。和我們常識(shí)理解的想法一致俊扭，GDP增加队橙，幸福指數(shù)會(huì)增加。
3）回歸系數(shù)的標(biāo)準(zhǔn)誤差stand error：0.061萨惑，即β的估計(jì)的標(biāo)準(zhǔn)差捐康。通過每年的GDP數(shù)據(jù)，可以計(jì)算得到回歸系統(tǒng)的標(biāo)準(zhǔn)誤差庸蔼〗庾埽回歸系數(shù)標(biāo)準(zhǔn)誤差，是量度結(jié)果精密度的指標(biāo)姐仅。這里計(jì)算得出的標(biāo)準(zhǔn)誤差為0.061花枫，數(shù)值很小刻盐，說明精確度還是不錯(cuò)的。
（ 補(bǔ)充：標(biāo)準(zhǔn)差是表示個(gè)體間變異大小的指標(biāo),反映了整個(gè)樣本對(duì)樣本平均數(shù)的離散程度,是數(shù)據(jù)精密度的衡量指標(biāo);而標(biāo)準(zhǔn)誤反映樣本平均數(shù)對(duì)總體平均數(shù)的變異程度,從而反映抽樣誤差的大小 ,是量度結(jié)果精密度的指標(biāo)劳翰。）
4）p-值為0%敦锌。根據(jù)簡單線性回歸顯著性的t檢驗(yàn)，原假設(shè)happiness 同gdp_per_capita之間不存在線性關(guān)系佳簸，β為0乙墙。而現(xiàn)在p值為0%，小于顯著性水平0.05生均。所以拒絕原假設(shè)听想，β顯著不等于0。我們足以斷定马胧，happiness 同gdp_per_capita之間存在一個(gè)顯著的關(guān)系汉买。
5）β的63%的置信區(qū)間：2.076 ~ 2.315。我們有63%的信心佩脊，回歸系數(shù)β將落在置信區(qū)間 [2.076蛙粘，2.315]中。換個(gè)角度來講邻吞，簡單線性回歸顯著性的t檢驗(yàn)组题，假設(shè)β為0，而β=0并沒有包含在上述置信區(qū)間內(nèi)抱冷，所以我們可以拒絕原假設(shè)崔列，斷定happiness 同gdp_per_capita之間存在一個(gè)顯著的關(guān)系。
3旺遮、多元線性回歸
回歸模型
我們知道僅僅考慮GDP是不能完全解釋幸福指數(shù)的赵讯。為了更加準(zhǔn)確地分析影響幸福指數(shù)因素，應(yīng)該引入一些不同的自變量來分析耿眉，看看哪些自變量的組合更優(yōu)地滿足OLS最小二乘假定边翼。這里利用六個(gè)自變量來測試一下。
from statsmodels.formula.api import ols
# 建立多元線性回歸模型
lm_m = ols(formula='happiness ~ gdp_per_capita + healthy_life_expectancy + freedom_to_life_choise + corruption_perceptions + generosity',
data=df_model).fit()
lm_m.summary()
1） 修正判定系統(tǒng)Adj.R-squared從之前的63.3%%上升到74.3%鸣剪。房屋價(jià)格指數(shù)變異性的74.3%能被其與多個(gè)自變量之間的線性關(guān)系解釋组底。
2）回歸系數(shù): 控制其他變量不變的情況下，GDP指數(shù)每增加一個(gè)單位筐骇，幸福指數(shù)增加1.32個(gè)單位瑞筐，健康預(yù)期壽命指數(shù)每增加一個(gè)單位庆械，幸福指數(shù)增加1.21個(gè)單位焕窝。
3）變量重要性排序?yàn)椋篻dp_per_capita棒呛、freedom_to_life_choise、healthy_life_expectancy告唆、corruption_perceptions棺弊、generosity
5）新引入的自變量的p值都小于顯著性水平0.05晶密，說明這些自變量同幸福指數(shù)是有顯著關(guān)系的∧Ｋ考慮到判定系數(shù)的增加稻艰，更加說明多元線性回歸模型在這里是優(yōu)于簡單一元線性回歸的。
六侈净、比較預(yù)測值和真實(shí)值的分布
y_pred = lm_m.predict(df_model[:])
df_pred = pd.concat([df_model['happiness'], y_pred], axis=1)
df_pred.columns = ['y_true', 'y_pred']
# 散點(diǎn)圖
fig = px.scatter(df_pred, x='y_true', y='y_pred', trendline='ols')
fig.update_layout(title='Resid of OLS Regression')
plot(fig, filename='./html/預(yù)測值和真實(shí)值分布圖.html')
七连锯、模型殘差分布圖
殘差: y預(yù)測值與y之間的差值
fig = px.histogram(x=lm_m.resid)
fig.update_layout(title='Resid of OLS Regression')
plot(fig, filename='./html/多元線性回歸殘差分布圖.html')
希望本文的內(nèi)容對(duì)大家的學(xué)習(xí)或者工作能帶來一定的幫助，每天進(jìn)步一點(diǎn)點(diǎn)用狱，加油~
