多元统计分析

# -*- coding: utf-8 -*-
"""
Created on Fri Apr 29 19:42:55 2022

@author: 25526
"""

#多元线性回归
#1.导入模块
#住房房价预测
import random
from pandas.plotting import scatter_matrix
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.model_selection import train_test_split
from sklearn.linear_model import Ridge,RidgeCV
from sklearn.linear_model import Lasso,LassoCV
from sklearn.metrics import mean_squared_error
import numpy as np
#from sklearn.model_selection import StratifiedKFold
###################任务一
##利用指定的12个自变量与因变量 comedy 创建散布图矩阵，主要目的查看各自变量与因变量之间的相关性。
#读取数据文件
pd_data1=pd.read_csv(r"./BostonHousing2.csv",header=1)
dff=pd_data1[['crim','zn','indus','nox','rm','age','dis','rad','tax','ptratio','b','lstat','cmedv']]
#plt.figure()#自变量
#scatter_matrix(frame,alpha=0.5,figsize=None,marker='.')
#frame,pandas dataframe对象
#alpha,图像透明度，一般取（0，1]
#figsize,以英寸为的单位的图像大小，一般以元组（width,heigth)形式设置
#marker,Matplotlib可用的标记类型，如'.',',','o'等
scatter_matrix(dff,alpha=0.5,figsize=(12,8))
#plt.grid()#生成网格
#plt.savefig('cmedy')
plt.show()
########################任务二
############回归分析
#随机地将当前数据集按照 3:1 的容量比例划分为训练集（用于建立模型）和测试集(用于检测模型的预测精度)，
#重复此步骤十次，并将得到十次结果制作如下的折线图，其中横坐标为次数，纵坐标为对应次数的可决系数
###交叉验证：测试集与训练集的划分
#X_train,X_test,Y_train,Y_test=train_test_split(X,Y,test_size=None,train_size=None,random_state=None,shuffle=True,stratify=None)
#X_train,X_test,Y_train,Y_test划分出的训练集/测试集 数据/标签
#X :待划分的像本特征集合（array-like）, Y : 待划分的样本标签
#test_size默认为None，若在0-1之间，为测试集样本数目与原始样本数目之比；若为整数，则是测试集样本的数目
#train_size默认为None,
#random_state随机数种子
X=pd_data1[['crim','zn','indus','nox','rm','age','dis','rad','tax','ptratio','b','lstat']]
Y=pd_data1[['cmedv']]
l=[]
for i in range(10):
    X_train,X_test,Y_train,Y_test=train_test_split(X,Y,train_size=0.75)
#2.创建模型对象
##最小二乘线性回归
#model=LinearRegression(fit_intercept=True,normalize=False,copy_X=True,n_jobs=1)
#fit_intercept默认为True，是否计算模型截距；为False时，则数据中心化处理
#normalize默认False,是否中心化，或使用sklearn.preprocessing.StandarScaler()
#copy_X=True,n_jobs默认为1，表示使用CPU的1个数。当-1时，代表使用全部CPU
    #model=LinearRegression()
    model = Ridge(alpha=0.5,fit_intercept=True)
    model = RidgeCV(alphas=[0.1, 1.0, 10.0]) # 通过 RidgeCV 使用交叉验证获取最佳参数值  
    model.fit(X_train,Y_train)
    a=model.intercept_ #截距
    b=model.coef_ #回归系数
#3.调用模型对象方法和属性实现
#a.模型训练    #train训练集，test测试集
    model.fit(X_train,Y_train)
#b.模型预测
    y_train_pred=model.predict(X_test)
#c.模型度量
    score=model.score(X_test,Y_test)#可决系数
    l.append(score)
####matplotlib绘制折线图
#1).获取自己所需点的坐标：横坐标为次数，纵坐标为对应次数的可决系数
x = ["第1次","第2次","第3次","第4次","第5次","第6次","第7次","第8次","第9次","第10次"]
y = [round(i,2) for i in l]
#图像上每个点上的数值：ha参数控制水平对齐方式, va控制垂直对齐方式
for a,b in zip(x,y):
    plt.text(a,b,b,ha='center',va='bottom')
#2).连线绘图,plot():绘制折线图的函数
#color:表示折线的颜色，label：折线的标签，
#linewidth：折线设置宽度（粗细），marker=“o”：可以使得折线中每个点被标记出来
plt.plot(x,y,linewidth=1, color='g', markerfacecolor='blue', marker='o')
#3).折线图坐标轴设置：调用label函数
plt.xlabel(u"次数")
plt.ylabel(u"可决系数")
#4).图像标题设置：
plt.title(u"随机划分下的模型可决系数")
#背景网格设置：plt.grid()
#折线说明框设置：plt.legend()
plt.ylim(0, 1) # 限定纵轴的范围
plt.show()

###########任务三##########
#(3)最优回归方程的选择:从12个自变量中随机的抽取n(其中n=2,…..12)个自变量，
#并利用十折交叉验证计算所建模型的可决系数，
#依据以上 12 个模型的可决系数大小确定哪一个模型的预测精度较高。（并不一定使用全部自变量的模型精度最好）
print("任务三")
print("下面对变量进行随机抽取：")
X_train,X_test,Y_train,Y_test = train_test_split(X,Y,train_size=0.75)
for p in range(10):
    ans=random.randint(1, 12)
    df3 = X_train.sample(n=ans,axis=1)
    l=list(df3.columns)
    df4=pd.DataFrame(X_test,columns=l)
    model = Ridge(alpha=0.5,fit_intercept=True)
    model = RidgeCV(alphas=[0.01,0.1,0.2, 0.5, 1],cv=10)
    model.fit(df3,Y_train)
    ridge_best_alpha = model.alpha_ #得到最佳lambda值
    a=model.intercept_
    b=model.coef_
    yy_train_pred =model.predict(df3)
    score = model.score(df4,Y_test)
    print('第{}轮,随机抽取{}个变量'.format(p+1,ans))
    print('岭回归关键正则参数={},可决系数{}'.format(ridge_best_alpha,round(score,2)))

################任务四##############
#（4）岭回归、Lasso 回归模型中关键正则参数𝛌的选择：在给定参数𝛌的 0.01，0.1，0.2, 0.5, 1 这五个可能性取值的条件下，利用十折交叉验证和可决系数确定两个模型的各自最优参数𝛌。
print("任务四")
print("岭回归回归模型中关键正则参数𝛌的选择:")
X_train,X_test,Y_train,Y_test = train_test_split(X,Y,train_size=0.75)
model = Ridge(alpha=0.5,fit_intercept=True)
model = RidgeCV(alphas=[0.01,0.1,0.2, 0.5, 1],cv=10)
model.fit(X_train,Y_train)
ridge_best_alpha = model.alpha_ #得到最佳lambda值
print('岭回归关键正则参数𝛌={}'.format(ridge_best_alpha))
ridge = Ridge(alpha = ridge_best_alpha,normalize = True)
ridge.fit(X_train,Y_train)
ridge_predict = ridge.predict(X_test)
mse = np.sqrt(mean_squared_error(Y_test,ridge_predict))
score=model.score(X_test, Y_test)
print('在最优参数的条件下,损失值{},可决系数:{}'.format(round(mse,2),round(score,2)))
print("LASSO回归模型中关键正则参数𝛌的选择:")
X_train,X_test,Y_train,Y_test = train_test_split(X,Y,train_size=0.75)
lasso_cv = LassoCV(alphas = [0.01,0.1,0.2, 0.5, 1], cv = 10)
lasso_cv.fit(X_train,pd.DataFrame(Y_train).values.ravel())  #本来需要一维向量，但是输入为列向量，所以找到出错的位置，应用ravel()函数即可。
lasso_best_alpha = lasso_cv.alpha_  # 取出最佳的lambda值
print('LASSO回归关键正则参数𝛌={}'.format(lasso_best_alpha))
lasso = Lasso(alpha = lasso_best_alpha, normalize=True)
lasso.fit(X_train,pd.DataFrame(Y_train).values.ravel())
lasso_predict = lasso.predict(X_test) #预测
MSE = np.sqrt(mean_squared_error(Y_test,lasso_predict))
ss=model.score(X_test, Y_test)
print('在最优参数𝛌的条件下,损失值{},可决系数:{}'.format(round(MSE,2),round(ss,2)))
###############任务五############
print("任务五")
ll=lasso.coef_
print("各变量回归系数如下:")
print(ll)
ans=0
for i in ll:
    if abs(i)<=0.05:
        ans+=1
tmp=round(ans/12,2)
tmp=int(tmp*100)
print('回归系数为零的个数与全部自变量个数比例:{}%.'.format(tmp))

#####Fisher判别分析##########

# -*- coding: utf-8 -*-
"""
Created on Wed May  4 06:37:38 2022

@author: 25526
"""
#############Fisher线性判别分析################
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
from sklearn.model_selection import train_test_split
from sklearn import discriminant_analysis 
from sklearn.model_selection import KFold
#（1）LDA 建模
#读取数据集合 pendigits.csv 文件(V17 为类标签)，并将其随机按照 1:5 的比例划分为训练集和测试集，估计模型的分类准确率。
pd_data=pd.read_csv("pendigits.csv")
X=pd_data.drop('V17',axis=1)
Y=pd_data[['V17']]
X_train,X_test,Y_train,Y_test=train_test_split(X,Y,train_size=0.84)
lda=discriminant_analysis.LinearDiscriminantAnalysis()
lda.fit(X_train,pd.DataFrame(Y_train).values.ravel())#消除警告
a=lda.coef_ #权向量值
b=lda.intercept_ #截距
lda.predict(X_train) 
#lda.score(X_test,Y_test)  #模型精度
print('Score: %.2f' % lda.score(X_test, Y_test))

#2.高维数据建模
#利用十折交叉验证和 Fisher 线性判别准则对数据集 LDA-sparse_data 进行建模，
#并观察指定模型中相关输入参数 shrinkage 的取值与不指定该参数取值时模型前后之间的差异。
pd_data1=pd.read_csv("LDA-sparse_data.csv")
#drop([,axis=0,inplace=False]):drop()删除行和列，默认删除某一行；若axis=1，删除某一列；inplace默认为False，保持原数据不变，True，改变
X1=pd_data1.drop('label',axis=1) 
Y1=pd_data1[['label']]
X1_train,X1_test,Y1_train,Y1_test=train_test_split(X1,Y1,train_size=0.8)
#np.linspace(start, stop, num=50, endpoint=True, retstep=False, dtype=None, axis=0)
#start:返回样本数据开始点,stop:返回样本数据结束点,num:生成的样本数据量，默认为50,endpoint：True则包含stop；False则不包含stop,retstep：如果为True则结果会给出数据间隔,dtype：输出数组类型,axis：0(默认)或-1
shrinkages=np.linspace(0.0,1.0,num=20)  #用来创建等差数列
scores=[]
for sh in shrinkages:
    #solver：逆矩阵求解方式，‘svd‘(默认), ’lsqr‘, ’eigen‘(数据稀疏时使用)
    #Shrinkage：收缩系数，当数据稀疏时配合’lsqr‘, ’eigen‘使用：’auto‘; float, between 0 and 1
    # n_components：投影后数据的维数，默认min(n_classes - 1, n_features)
    lda1=discriminant_analysis.LinearDiscriminantAnalysis(solver='lsqr',shrinkage=sh)
    #pd.DataFrame()创建二维表,df_obj.values #查看数据值
    #np.ravel()将多维数组降为一维,返回的是视图,会影响原始矩阵。
    lda1.fit(X1_train,pd.DataFrame(Y1_train).values.ravel()) #训练模型
    lda1.predict(X1_train)
    scores.append(lda1.score(X1_test,Y1_test))
print("指定参数,shrinkage从[0,1]取20个,Score分别为：{}".format(scores))
lda2=discriminant_analysis.LinearDiscriminantAnalysis()
lda2.fit(X1_train,pd.DataFrame(Y1_train).values.ravel())
lda2.predict(X1_train)
print("不指定参数，score取值：{}".format(lda2.score(X1_test,Y1_test)))

#任务三
#3. 模型的适用性
#（1）读取数据集 banana.dat，并在二维坐标系中绘制该数据集的散点图，其中两类数据分别使用两种不同颜色的点表示；
#（2）利用十折交叉验证和 Fisher 线性判别法对该数据集进行建模，并观察分类准确率是否较高
data2=pd.read_csv(r'./banana.dat')
X2=data2.drop(' Class',axis=1)
Y2=data2[[' Class']]
#作用：将当前的数据划分按照既定的比例划分为若干个训练集和测试集的组合，并返回模型对象以供使用
#kfolds = KFold(n_splits=5, shuffle=False, random_state=None)
#n_splits :折的数量;shuffle ：是否打乱顺序;random_state：随机数种
folder=KFold(n_splits=10,shuffle=False,random_state=None)
#Kfolds.split(X,y): 对数据集X进行切分;
#输入：数据array的X,y;输出：返回测试集与训练集的K个组合的索引值
avg_folder=0
for train,test in folder.split(X2,Y2):
    X2_train=X2.iloc[train,:]
    X2_test=X2.iloc[test,:]
    Y2_train=Y2.iloc[train,:]
    Y2_test=Y2.iloc[test,:]
    lda3=discriminant_analysis.LinearDiscriminantAnalysis()
    lda3.fit(X2_train,pd.DataFrame(Y2_train).values.ravel())
    lda3.predict
    avg_folder+=lda3.score(X2_test,Y2_test)
avg_folder=avg_folder/10
print("用十折交叉验证和 Fisher 线性判别法,准确率为{}".format(avg_folder))

print("***************************散点图如下*****************")
#sns.scatterplot(x,y,hue,size,style)
#x,y指定x轴和y轴位置上的变量；
#hue：对将生成具有不同颜色的点的变量进行分组；
#size：将生成不同大小的点的变量进行分组；
#style：将生成具有不同标记的点的变量进行分组
sns.scatterplot(data2['At1'],data2[' At2'],
                hue=data2[' Class'],size=data2[' Class'],
                style=data2[' Class'])  #右图，加上hue参数

########主成分分析###########

# -*- coding: utf-8 -*-
"""
Created on Sat May  7 10:26:21 2022

@author: 25526
"""
import pandas as pd
import numpy as np
import sklearn.decomposition as dp
import matplotlib.pyplot as plt

data=pd.read_csv(r'./secom.data',delimiter=' ') #delimiter 自定义分隔符
num=data.shape[1] #shape[0]的矩阵的行数；shape[1]列数
#######案例一、半导体数据降维##########
#1.数据预处理:本数据集存在缺失数据，将每列的缺失值补全为该列的所有非缺失值的均值。
for column in list(data.columns[data.isnull().sum()>0]):
    mean_val=data[column].mean()   #mean()求取每列特征均值
    #fillna()填充缺失值，
#    #若inplace为True，直接修改原对象，若为False，默认值，创建一个副本，修改副本，原对象不变
    data[column].fillna(mean_val,inplace=True)  
#2.利用 sklearn 模块相关方法进行主成分分析
#（1）画出崖底碎石图，观察此图看是否个主成分的贡献率的差异情况；
#data1=data.T #对数据行列进行转置
#ov_mat = (data1.T - np.mean(data1 , axis = 1)).T#  x-mean(x) 以备计算协方差矩阵使用
#eig_val, eig_vec = np.linalg.eigh(np.cov(data1)) #计算协方差矩阵的特征值和特征向量
#p = np.size(eig_vec, axis =1)
#idx = np.argsort(eig_val) #按照特征值从小到大顺序排序
#idx = idx[::-1]#按照特征值从大到小顺序排序
#eig_vec = eig_vec[:,idx] #特征向量，即计算主成分的权向量（所有的）
#eig_val = eig_val[idx]      #特征值，即所有主成分的方差
#ll=eig_val/np.sum(eig_val)
#plt.rcParams['font.sans-serif']=['SimHei']
#plt.rcParams['axes.unicode_minus'] = False
#plt.plot(ll)
#plt.xlabel("主成分的个数")
#plt.ylabel("主成分贡献率")
#plt.show()
#（2）选择不同的阈值 0.7，0.8，0.9，筛选出不同个数的主成分；
#（3）计算在不同阈值条件下属性的压缩比，即主成分个数/所有原始数据的属性个数
n_com=[0.7,0.8,0.9]
for n_component in n_com:
    pca=dp.PCA(n_components=n_component)
    reduced_x=pca.fit_transform(data) 
    reduced_x= np.dot(data,pca.components_.T)
    x=reduced_x.shape[1]
    ratio=round(x/num,5)
    print("在阈值{}下,主成分个数{}个,压缩比为{}".format(n_component,x,ratio))
    

#######案例二、图像特征抽取##############
#1. 数据采集与处理:
#自拍一张照片，然后调用 PCA-fig.py 文件中的相关函数，并制定主成分的个数 1,2,3…100，
#观察一下在不同主成分个数的条件下，图片的变化。
#2. 程序修改:设定主成分累积贡献率的阈值实现主成分个数的自适应确定