第4章 探索性數據分析（多因子與復合分析）

4.1、多因子

4.1.1 假設檢驗與方差檢驗

　　假設檢驗適用於（數據樣本較小時）

　　方差檢驗適用於（數據樣本較大時）

import numpy as np
import scipy.stats as ss
#生成一20個數的標准正態分布
norm_dist = ss.norm.rvs(size=20)
#檢測norm_dist是否是正態分布，使用的方法是基於峰度和偏度的
print(ss.normaltest(norm_dist))

#結果：NormaltestResult(statistic=0.2025598777545946, pvalue=0.9036800223028876)
     #第一個是統計值，第二個值是p值

（1）P分布檢驗常用於比較兩種樣本是否一致（例如：臨床醫療上葯物是否有效）；

（2）獨立分布t檢驗用於檢測兩組值的均值是都有比較大的差異性

print(ss.ttest_ind(ss.norm.rvs(size=10),ss.norm.rvs(size=20)))
#結果：Ttest_indResult(statistic=-0.575484958550556, pvalue=0.5695598474341583)

　　由於p值大於0.05（假定），可以接受該假設

（3）卡方檢驗常常用於確定兩因素件是否有比較強的聯系

import scipy.stats as ss
print(ss.chi2_contingency([[15,95],[85,5]]))
#結果：(126.08080808080808, 2.9521414005078985e-29, 1, array([[55., 55.],[45., 45.]]))
     #第一個值為卡方值，第二個值為P值，第三個值為自由度，第四個為與原數據數組同維度的對應理論值

　由於p值小於0.05（假定），拒絕該假設

（4）F檢驗常用於方差分析

　　　檢測統計量F，做假設檢驗【F滿足自由度（m-1，n-m）的F分布】

SST：總變差平方和

SSM：平均平方和（組間平方和）

SSE：殘差平方和（組內平方和）

print(ss.f_oneway([49,50,39,40,43],[28,32,30,26,34,],[38,40,45,42,48]))
#結果：F_onewayResult(statistic=17.619417475728156, pvalue=0.0002687153079821641)

（5）使用qq圖來對比一個分布是否屬於一個已知的分布

  　　　方法：看散點圖的是否在x，y軸的角平分線上

from statsmodels.graphics.api import qqplot
import matplotlib.pyplot as plt
qqplot(ss.norm.rvs(size=100))
plt.show()

 

4.1.2 相關系數

范圍在[-1,1]----皮爾遜、斯皮爾曼（只與名次差有關，與具體數值相關性不強，故適用於相對比較的情況下）

　　衡量兩組數據（或者兩組樣本）的分布趨勢、變化趨勢、一致趨勢程度的因子，分為正相關（正向同步），負相關（反向同步）等。

　　　　其中：分子是協方差，分母是x，y的標准差，兩者相乘相當於一個歸一化因子

 

　　　　其中：n表示每組數據的數量，d值得是兩組數據排名后的名次差

import pandas as pd
x = pd.Series([2.5,0.5,2.2,1.9,3.1,2.3,2,1,1.5,1.1])
y = pd.Series([2.4,0.7,2.9,2.2,3.0,2.7,1.6,1.1,1.6,0.9])
data = pd.DataFrame({'x':x,'y':y})
print('data:','\n',data)
print('Pearson相關系數：',x.corr(y))
print('Spearman相關系數：',x.corr(y,method='spearman'))
#Pearson相關系數： 0.9259292726922455
#Spearman相關系數： 0.9179373709568976

#DataFrame是針對於列進行相關性計算的，故需要對其進行轉置
df = pd.DataFrame(np.array([x,y]).T)
print(df.corr())
'''結果
          0         1
0  1.000000  0.925929
1  0.925929  1.000000
'''

4.1.3 回歸

線性回歸

　　線性回歸最常用的是解法：最小二乘法，其本質是最小化誤差平方

 

線性回歸效果判定的關鍵指標：

　　決定系數越接近於1說明效果越好；

　　DW的范圍在【0 , 4】，DW=2時表示殘差不相關，接近於4表示正相關，接近於0表示殘差負相關，好的回歸一般殘差是不相關的。

  其中決定系數中的左R²：一元決定系數的定義；

　　　　　　　　 右R²：多元決定系數的定義；

　　　　　　　　　k：參數的個數

　　　　　　　　　e：殘差（預測值與實際值的差）

https://scikit-learn.org/stable/modules/generated/sklearn.linear_model.LinearRegression.html#sklearn.linear_model.LinearRegression

import numpy as np
import scipy.stats as ss

x = np.arange(10).astype(np.float).reshape((10,1))
y =x*3+4+np.random.random((10,1))
print('實際值：','\n',y)
from sklearn.linear_model import LinearRegression
reg = LinearRegression()
res = reg.fit(x,y)
y_pre = reg.predict(x)
print('預測值：','\n',y_pre)
score = reg.score(x,y)
print('效果：',score)
a = reg.coef_
b = reg.intercept_
print('參數：',a,';截距：',b)
import seaborn as sns
import matplotlib.pyplot as plt
plt.figure(figsize=(10,8))
plt.plot(x,y)
plt.plot(x,y_pre)
plt.show()

　　

4.1.4 主成分分析（PCA）與奇異值分解（SVD）

主成分分析（PCA）步驟

·         求特征協方差矩陣

·         求協方差的特征值和特征向量

·         將特征值從大到小排序，選擇其中最大的k個

·         將樣本點投影到選取的特征向量上

https://scikit-learn.org/stable/modules/generated/sklearn.decomposition.PCA.html#examples-using-sklearn-decomposition-pca

自寫版本（按照PCA步驟）：

　　特別說明：np.linalg.eig(cov)之后的特征向量已按照對應的特征值從大至小進行排列

　　numpy.linalg模塊包含線性代數的函數。使用這個模塊，可以計算逆矩陣、求特征值、解線性方程組以及求解行列式等
    　　eigvals函數可以計算矩陣的特征值，而eig函數可以返回一個包含特征值和對應的特征向量的元組

 #-----特別說明：np.linalg.eig(cov)之后的特征向量已按照對應的特征值從大至小進行排列
import numpy as np
import pandas as pd
x=[2.5,0.5,2.2,1.9,3.1,2.3,2,1,1.5,1.1]
y=[2.4,0.7,2.9,2.2,3.0,2.7,1.6,1.1,1.6,0.9]
data = pd.DataFrame({'x':x,'y':y})
print('data:','\n',data)
cov = data.cov()
print('協方差矩陣:','\n',cov)
#求協方差矩陣的特征值和特征向量
#numpy.linalg模塊包含線性代數的函數。使用這個模塊，可以計算逆矩陣、求特征值、解線性方程組以及求解行列式等
    #eigvals函數可以計算矩陣的特征值，而eig函數可以返回一個包含特征值和對應的特征向量的元組
eigen = np.linalg.eig(cov)
print(eigen)
eigenvalues = eigen[0]
print('eigenvalues:','\n',eigen[0])
#選取比較大的特征值所代表的特征向量
eigenvectors = eigen[1]
print('選取的特征值:','\n',eigen[0][1])
print('選取的特征向量:','\n',eigen[1][0])
jian = data-np.mean(data)
pcavectors = np.dot(jian,eigen[1][0])
print('主成分分析后的向量','\n',pcavectors)

 自寫版本（使用函數按照PCA步驟）：

import numpy as np
import pandas as pd
def myPCA(data,n_components =1000000000):
    mean_vals = np.mean(data,axis=0)
    mid = data - mean_vals
    #此處需要針對列進行協方差計算
    cov_mat = np.cov(mid,rowvar=False)#rowvar默認為True，即對行求協方差
    from scipy import linalg
    eig_vals,eig_vects = linalg.eig(np.mat(cov_mat))
    #取出最大的特征值所對應的特征向量
    eig_val_index = np.argsort(eig_vals)#argsort()得到排序之后的索引
    eig_val_index = eig_val_index[:-(n_components+1):-1]
    eig_vects = eig_vects[:,eig_val_index]
    low_dim_mat = np.dot(mid,eig_vects)
    return low_dim_mat,eig_vals
x=[2.5,0.5,2.2,1.9,3.1,2.3,2,1,1.5,1.1]
y=[2.4,0.7,2.9,2.2,3.0,2.7,1.6,1.1,1.6,0.9]
data = pd.DataFrame([x,y]).T
print(myPCA(data,n_components=1))

接口版本：

from sklearn.decomposition import PCA
import numpy as np
import pandas as pd
x=[2.5,0.5,2.2,1.9,3.1,2.3,2,1,1.5,1.1]
y=[2.4,0.7,2.9,2.2,3.0,2.7,1.6,1.1,1.6,0.9]
data = pd.DataFrame([x,y]).T
#降維
lower_dim = PCA(n_components=1)  #降為1維
new = lower_dim.fit(data)
print(new)
　　#PCA(copy=True, iterated_power='auto', n_components=1, random_state=None,
　　#svd_solver='auto', tol=0.0, whiten=False)
　　#copy如果為False，則傳遞給fit的數據將被覆蓋並且運行fit（X）.t
important = lower_dim.explained_variance_ratio_
print(important)  #[0.96318131],說明得到了96.3%的信息量
newdata = lower_dim.transform(data)
print(newdata)    #轉化后的數值

注意：sklearn中所用到的PCA方法是奇異值分解（SVD），與上述PCA的步驟不一樣，所以兩者所得出的結論會有些許不同

4.2、復合分析

4.2.1交叉分析

　　除了縱向分析橫向分析外，可以使用交叉分析來分析數據屬性和屬性的關系

sns.heatmap()用於做熱圖-----兩個不同屬性間相互影響的情況！

import numpy as np
import pandas as pd
import scipy.stats as ss
import matplotlib.pyplot as plt
import seaborn as sns
df = pd.read_csv("HR.csv")
##查看各個部門的離職之間是否有明顯差異----使用獨立t檢驗方法
    ##得到各個部門的離職分布

#得到分組后的索引值
dp_indices = df.groupby(by = 'department').indices
# print(dp_indices)
sales_values = df["left"].iloc[dp_indices['sales']].values
technical_values = df["left"].iloc[dp_indices['technical']].values
print(ss.ttest_ind(sales_values,technical_values))

#求兩兩部門間的p值
dp_keys = list(dp_indices.keys())
print(dp_keys)
dp_t_mat = np.zeros([len(dp_keys),len(dp_keys)])
for i in range(len(dp_keys)):
    for j in range(len(dp_keys)):
        p_values = ss.ttest_ind(df['left'].iloc[dp_indices[dp_keys[i]]].values,
                                df['left'].iloc[dp_indices[dp_keys[j]]].values)[1]
        if p_values<0.05:
            dp_t_mat[i][j] = -1
        else:
            dp_t_mat[i][j] = p_values
sns.heatmap(dp_t_mat,xticklabels=dp_keys,yticklabels=dp_keys)
plt.show()
#只有為黑色的，兩個部門的才有顯著性差異

 

# 采用透視表
piv_tb = pd.pivot_table(df,values='left',index=['promotion_last_5years','salary'],
                        columns=['Work_accident'],aggfunc=np.mean)
print(piv_tb)
sns.heatmap(piv_tb,vmin=0,vmax=1,cmap=sns.color_palette("Reds",n_colors=256))#指定最大值為1，最小值為0
plt.show()

顏色越深，離職率越高

4.2.2分組與鑽取

鑽取是改變數據維度層次，變換分析粒度的過程，可分為向下和向上鑽取。
連續屬性在分組前需要進行離散化。
連續分組：

·         分割（一階差分）、拐點（二階差分）

·         聚類

·         不純度（Gini系數）

#求Gini系數
s1 = pd.Series(['X1','X1','X2','X2','X2'])
s2 = pd.Series(['Y1','Y1','Y1','Y2','Y2'])
#求平方和：
def getProbSS(s):
    if type(s)!= pd.core.series.Series:
        s = pd.Series(s)
    prt_ary = s.groupby(by=s).count().values/float(len(s))
    return 1- sum(prt_ary**2)
def getGini(s1,s2):
    d = dict()
    for i in list(range(len(s1))):
        d[s1[i]] = d.get(s1[i], []) + [s2[i]]
    return 1 - sum([getProbSS(d[k]) * len(d[k])/float(len(s1)) for k in d])
print("Gini:",getGini(s1,s2))
#結果：Gini: 0.7333333333333334

　

sns.barplot(x='salary',y='left',hue='department',data=df)
plt.show()

#拐點（二階差分）
sl_s = df['satisfaction_level']
sns.barplot(list(range(len(sl_s))),sl_s.sort_values())
plt.show()

4.2.3相關分析

#heatmap()會自動屏蔽離散化信息。例如，salary有low、medium、high

sns.heatmap(df.corr(),vmin=-1,vmax=1,cmap=sns.color_palette("RdBu",n_colors=128))
plt.show()

對離散屬性相關性的處理：

　　利用熵來進行離散屬性相關性的計算

        　　熵：用於衡量不確定性，條件熵不具有對稱性

 

s1 = pd.Series(['X1','X1','X2','X2','X2','X2'])
s2 = pd.Series(['Y1','Y1','Y1','Y2','Y2','Y2'])
def getEntropy(s):
    if type(s)!= pd.core.series.Series:
        s = pd.Series(s)
    prt_ary = s.groupby(by=s).count().values/float(len(s))
    return (-(np.log2(prt_ary))*prt_ary).sum()
#條件熵需要注意
def getCondEntropy(s1,s2):#s1條件下s2的熵
    d = dict()
    for i in list(range(len(s1))):
        d[s1[i]] = d.get(s1[i],[])+[s2[i]]
        #{'X1': ['Y1', 'Y1'], 'X2': ['Y1', 'Y2', 'Y2', 'Y2']}
    return sum([getEntropy(d[k]) * len(d[k])/float(len(s1)) for k in d])
def getEntropyGain(s1,s2):
    return (getEntropy(s2) - getCondEntropy(s1,s2))
def getEntropyGainRatio(s1,s2):
    return getEntropyGain(s1,s2)/getEntropy(s2)
import math
def getCoor(s1,s2):
    return getEntropyGain(s1,s2)/math.sqrt(getEntropy(s1)*getEntropy(s2))
print('熵:',getEntropy(s2))
print('條件熵:',getCondEntropy(s1,s2))
print('熵增益:',getEntropyGain(s1,s2))
print('熵增益率:',getEntropyGainRatio(s2,s1))
print('相關性:',getCoor(s1,s2))
'''
結果：
熵: 1.0
條件熵: 0.5408520829727552
熵增益: 0.4591479170272448
熵增益率: 0.5
相關性: 0.4791387674918639
'''

4.2.4因子分析

　　探索性因子分析：通過協方差矩陣、相關性矩陣等指標分析多元屬性變量的本質結構，並可以進行轉化、降維等操作，得到影響目標屬性的最主要因子。

　　驗證性因子分析：驗證一個因子與關注的屬性之間是否有一定的關聯，有何關聯。

 

from sklearn.decomposition import PCA
my_pca = PCA(n_components=7)
#擬合時不允許出現離散型數據，drop默認情況下從行進行刪除
low_mat = my_pca.fit_transform(df.drop(labels=['salary','department','left'],axis=1))
#通過PCA降維后信息的准確性
print("Ratio:",my_pca.explained_variance_ratio_)
sns.heatmap(pd.DataFrame(low_mat).corr(),vmin=-1,vmax=1,cmap = sns.color_palette("RdBu",n_colors=128))
plt.show()
#結果：Ratio: [9.98565340e-01 8.69246970e-04 4.73865973e-04 4.96932182e-05 2.43172315e-05 9.29496619e-06 8.24128218e-06]

通過PCA，將原來的特征空間變成了一個正交的特征空間

4.2.5聚類分析、回歸分析（暫略）