批量椒盐噪声 python_linear 、rigde、lasso（python实现）

柔光的暖阳◎ 2022-11-19 11:23 41阅读 0赞

![a6fca30027ebe230129bc2ed1109c7be.png][]

linear Regression、rigde、lasso的python代码实现。

包括了若干算法，如**Adam下降**求解线性回归、ridge，**迭代ridge**、**坐标下降**求解lasso等。

代码github自取(点个star支持下 ？)

https://github.com/vincen-github/Machine-Learning-Codegithub.com

# author:vincen
    # date:2020/4/10
    # encoding:utf-8
    import numpy as np
    import pandas as pd
    import matplotlib.pyplot as plt
    import seaborn as sns
    
    #Parent Class
    class LinearRegression(object):
        '''
        LinearRegression Object(parent class):
            attribute:
                c:噪声的影响强度,取值[0,1] 默认为0
                shape:生成数据集的形状,默认为None     dtype = [sample_num,feature_num]
                mean:x所服从正态分布的均值向量       默认为None
                cov:x所服从正态分布的协方差矩阵      默认为None
                noise_miu：噪声所服从的正态分布的均值
                noise_sigma:噪声所服从的正态分布的方差
                noise默认服从标准正态分布
            function:
                generate_data(): 按照实例化对象的初始属性生成数据集
                regression_cond():计算线性回归中 X 的条件数
                regression_correlation():可视化矩阵X的相关系数矩阵
                fit():最小二乘估计训练线性模型(直接利用(X^T*X)^(-1)*X*y计算)，return alpha和RSS
                score():通过训练好的model计算测试误差
                predict():利用已训练好的model预测测试样本
                batch_gradient_descent():批量梯度下降 每一个epoch用整个训练集计算梯度
                stochastic_gradient_descent():随机梯度下降 每一个epoch中运用每一个样本计算梯度
                mini_batch_gradient_descent():小批量梯度下降 每一个epoch中运用等大小的训练集的子集计算梯度
                AdamOptimizer():Adaptive Moment Estimation
                square_R():goodness of fit
                modified_square_R():modified_square_R
        '''
        def __init__(self,c = 0,shape = None,mean = None,cov = None,noise_miu = 0,noise_sigma = 1):
                self.c = c
                self.shape = shape
                self.mean = mean
                self.cov = cov
                self.noise_miu = noise_miu
                self.noise_sigma = noise_sigma
    
        def generate_data(self,beta,random_state = None):
            """
            parameters:
                input:beta
                    The parameters of data generated,the last value of beta is the bias of generative dataset.
                notation:该方法会自动生成偏置全1列
                return : the generative dataset
                    type: DataFrame
            """
            #检查beta的维数是否是特征数+1(+1是因为偏置的存在)
            if len(beta) != (self.shape)[1]+1:
                raise Exception("The dimensions of beta is not fit the shape of dataset.")
            # 检查协方差矩阵和均值向量的维数，若否，抛出异常
            elif np.asarray(self.cov).shape[0] != self.shape[1]:
                raise Exception("The dimensions of the covariance matrix and the shape of dataset are not fitting.")
            elif len(self.mean) != self.shape[1]:
                raise Exception("The dimensions of the mean vector and the shape of dataset are not fitting")
            # 检查协方差矩阵是否为半正定阵，若否，抛出异常
            elif ~(np.all(np.linalg.eigvals(self.cov) >= 0 )):
                raise Exception('Covariance matrix is not positive define')
            else:
                # print(len(beta))
                np.random.seed(random_state)
                #利用标准正态分布生成X
                X = np.random.multivariate_normal(mean = self.mean, cov = self.cov, size = self.shape[0])
                # 将全1的列向量加入到X的最后1列
                X = np.hstack((X,np.ones(shape = (self.shape[0],1))))
                # print("the size of x:",self.X.shape)
                #加入噪声生成y
                y = np.matmul(X,np.asarray(beta).reshape([-1,1])) + self.c*np.random.normal(self.noise_miu,self.noise_sigma,size = (self.shape[0],1))
                # print("the size of y:",self.y.shape)
                #合并
                dataset = pd.DataFrame(np.hstack([X,y]))
                # return 
                return dataset
            
        def cond(self,A):
            '''
            parameters:
                input:
                    A: the coefficient matrix of Ax = b
            return:
                the condition number of A
            '''
            #  Do SVD decomposition on X to get singular values
            A = np.asarray(A)
            U,Sigma,_ = np.linalg.svd(A)
            # The condition number in L_2 space is equal to the largest singular value divided by the smallest singular value
            return max(Sigma)/min(Sigma)
    
        def correlation(self,A):
            '''
            parameters:
                input:
                    A: the coefficient matrix of Ax = b
                    dtype:DataFrame
            output:
                the figure of correlation matrix
            return:
                the correlation matrix of X
            '''
            # calcuate the correlation matrix of A,the last column,all values are 1,have been dropped by me.
            self.corrmat = A.corr()
            plt.subplots(dpi = 110)
            sns.heatmap(self.corrmat,vmin = -1,vmax = 1,square = True,cmap = 'RdYlBu_r',annot = True,fmt='.2f',annot_kws = {
         "size":8})
            plt.title('$Correlation quad Matrix$')
            plt.xticks(list(range(A.shape[1])),A.columns)
            plt.show()
            return self.corrmat
            
    
        def fit(self,dataset):
            """
            parameters:
                input:dataset
                    type:DataFrame
                        the last columns is label
                return : the parameters of linear model,traing_error
                    type: ndarray
            """
            X = np.asarray(dataset.iloc[:,:-1]) 
            y = np.asarray(dataset.iloc[:,-1]).reshape(-1,1)
            self.alpha = np.matmul(np.linalg.inv(np.matmul(X.T,X)),np.matmul(X.T,y))
            y_predict = np.matmul(X,self.alpha)
            self.training_error = ((y - y_predict)**2).sum()
            return self.alpha,self.training_error
    
        def score(self,testing_data):
            """
            parameters:
                input:testing_data
                    type:DataFrame
                        the last columns is label
                return : testing_error
                    type: ndarray
            """
            X_test = np.asarray(testing_data.iloc[:,:-1])
            print('1',X_test.shape)
            y_test = np.asarray(testing_data.iloc[:,-1]).reshape(-1,1)
            print('2',y_test.shape)
            y_predict = np.dot(X_test,self.alpha)
            print('3',y_predict.shape)
            print('4',(y_test.reshape(y_predict.shape) - y_predict).shape)
            self.testing_error = ((y_test.reshape(y_predict.shape) - y_predict)**2).sum()
            return self.testing_error
    
        def predict(self,dataset):
            """
            parameters:
                input:
                    dataset  d
                        type:DataFrame
                            dataset is not include labels
                return : the predict result
                    type: ndarray
            """
            return np.matmul(np.asarray(dataset),self.alpha)
    
        def batch_gradient_descent(self,dataset,max_epoch,eta = 0.0001,show_result = False):
            """
            parameters:
                input:
                    max_epoch    type:int
                    eta    type:float   learning rate Default:0.001
                    dataset    type:DataFrame    the last columns is label
                    show_result    type:bool    Default:False
                output:
                    if show_result == True:
                        output the figure of the process of training
                return : the parameters of linear model,traing_error
                    type: ndarray
            """
            self.eta = eta
            X = np.asarray(dataset.iloc[:,:-1])
            y = np.asarray(dataset.iloc[:,-1]).reshape(-1,1)
            self.error_list = []
            # sample the initial params from the unifrom distribution 
            self.alpha = np.random.uniform(0,1,size = [dataset.shape[1] - 1,1])
            # print("##########################",self.alpha.shape)
            for epoch in range(max_epoch):
                error = self.score(dataset)
                self.error_list.append(error)
                self.alpha -= self.eta*np.matmul(-X.T, y - np.matmul(X,self.alpha))
            if show_result == True:
                sns.set()
                plt.figure(dpi = 110)
                plt.plot(self.error_list,c = 'b',label = '$loss quad curve$')
                plt.title('$Batch quad Gradient quad Descent quad error$')
                plt.xlabel('$epoch$')
                plt.ylabel('$error$')
                plt.legend()
                plt.show()
            return self.alpha,self.error_list
        
        def AdamOptimizer(self,dataset,max_epoch,eta = 0.001,show_result = False):
            """
            parameters:
                input:
                    epsilon     type:float      the stop condition      Default:1e-6 
                    eta    type:float   learning rate   Default:0.001
                    dataset    type:DataFrame    the last columns is label
                    show_result    type:bool    Default:False
                output:
                    if show_result == True:
                        output the figure of the process of training
                return : the parameters of linear model,traing_error
                    type: ndarray
            """
            beta1 = 0.9
            beta2 = 0.999
            X = np.asarray(dataset.iloc[:,:-1])
            y = np.asarray(dataset.iloc[:,-1]).reshape(-1,1)
            self.error_list = []
            # sample the initial params from the unifrom distribution 
            self.alpha = np.random.uniform(0,1,size = [dataset.shape[1] - 1,1])
            # print("##########################",self.alpha.shape)
            error = self.score(dataset)
            m = np.zeros(shape = [X.shape[1],1])
            v = np.zeros(shape = [X.shape[1],1])
            t = 1
            while t < max_epoch:
                error = self.score(dataset)
                self.error_list.append(error)
                grad = np.matmul(-X.T, y - np.matmul(X,self.alpha))
                m = beta1*m + (1 - beta1)*grad
                v = beta2*v + (1 - beta2)*np.power(grad,2)
                m_hat = m/(1 - beta1**t)
                v_hat = v/(1 - beta2**t)
                self.alpha -= eta/(np.sqrt(v_hat)+1e-8)*m_hat
                t += 1
            
            if show_result == True:
                sns.set()
                plt.figure(dpi = 110)
                plt.plot(self.error_list,c = 'b',label = '$loss quad curve$')
                plt.title('$Adam quad Optimizer quad error$')
                plt.xlabel('$epoch$')
                plt.ylabel('$error$')
                plt.legend()
                plt.show()
            return self.alpha,self.error_list
    
        def square_R(self,dataset):
            """
            parameters:
                dataset:
                    the dataset which we want to compute  R^2 above it
                    dtype: DataFrame            the last column is label
                return:
                    R_square
            """
            SSE = self.score(dataset)
            SSR = np.power((self.predict(dataset.iloc[:,:-1]) - dataset.iloc[:,-1].mean()),2).sum()
            SST = SSE + SSR
            #注意这里已经把全1列加入到dataset中，所以不需要减1
            R_square = 1 - SSE/SST
            return R_square
        
        def modified_square_R(self,dataset):
            """
            parameters:
                dataset:
                    the dataset which we want to compute modified R^2 above it
                    dtype: DataFrame            the last column is label
                return:
                    modified R_square
            """
            SSE = self.score(dataset)
            SSR = np.power((self.predict(dataset.iloc[:,:-1]) - dataset.iloc[:,-1].mean()),2).sum()
            SST = SSE + SSR
            #注意这里已经把全1列加入到dataset中，所以不需要减1
            modified_R_square = 1 - (SSE/(dataset.shape[0] - dataset.iloc[:,:-1].shape[1] - 1))/(SST/(dataset.shape[0] - 1))
            return modified_R_square
    
    
    # Child Class RidgeRegression
    class RidgeRegression(LinearRegression):
        '''
        RidgeRegression Object(child class):
            attribute:
                c:噪声的影响强度,取值[0,1] 默认为0
                shape:生成数据集的形状,默认为None     dtype = [sample_num,feature_num]
                mean:x所服从正态分布的均值向量       默认为None
                cov:x所服从正态分布的协方差矩阵      默认为None
                noise_miu：噪声所服从的正态分布的均值
                noise_sigma:噪声所服从的正态分布的方差
                noise默认服从标准正态分布
            function:
                generate_data(): 按照实例化对象的初始属性生成数据集
                regression_cond():计算Ridege中 X 的条件数
                regression_correlation():可视化矩阵X的相关系数矩阵
                fit():最小二乘估计训练线性模型((直接利用(X^T*X)^(-1)+CI)*X*y计算)，return alpha和RSS
                score():通过训练好的model计算测试误差
                predict():利用已训练好的model预测测试样本
                batch_gradient_descent():批量梯度下降 每一个epoch用整个训练集计算梯度
                stochastic_gradient_descent():随机梯度下降 每一个epoch中运用每一个样本计算梯度
                mini_batch_gradient_descent():小批量梯度下降 每一个epoch中运用等大小的训练集的子集计算梯度
                AdamOptimizer():Adaptive Moment Estimation
                square_R():goodness of fit
                modified_square_R():modified_square_R
                trace():Ridge coefficients as a function of the regularization
        '''
        def __init__(self,c = 0,shape = None,mean = None,cov = None,noise_miu = 0,noise_sigma = 1,C = 1):
            super(RidgeRegression,self).__init__(c,shape,mean,cov,noise_miu,noise_sigma)
            self.C = C
        
        #重写父类LinearRegression的fit方法,变为RidgeRegression
        def fit(self,dataset):
            """
            parameters:
                input:dataset
                    type:DataFrame
                        the last columns is label
                return : the parameters of linear model,traing_error
                    type: ndarray
            """
            X = np.asarray(dataset.iloc[:,:-1])
            y = np.asarray(dataset.iloc[:,-1]).reshape(-1,1)
            self.alpha = np.matmul(np.linalg.inv(np.matmul(X.T,X)+self.C*np.diag(np.ones(X.shape[1]))),np.matmul(X.T,y))
            y_predict = np.matmul(X,self.alpha)
            self.training_error = ((y.reshape(y_predict.shape) - y_predict)**2).sum()
            return self.alpha,self.training_error
    
        #重写父类LinearRegression的batch_gradient_descent()方法
        def batch_gradient_descent(self,dataset,max_epoch,eta = 0.0001,show_result = False):
            """
            parameters:
                input:
                    max_epoch    type:int
                    eta    type:float   learning rate Default:0.0001
                    dataset    type:DataFrame    the last columns is label
                    show_result    type:bool    Default:False
                output:
                    if show_result == True:
                        output the figure of the process of training
                return : the parameters of linear model,traing_error
                    type: ndarray
            """
            self.eta = eta
            X = np.asarray(dataset.iloc[:,:-1])
            y = np.asarray(dataset.iloc[:,-1]).reshape(-1,1)
            self.error_list = []
            # sample the initial params from the unifrom distribution 
            self.alpha = np.random.uniform(0,1,size = [dataset.shape[1] - 1,1])
            # print("##########################",self.alpha.shape)
            for epoch in range(max_epoch):
                error = self.score(dataset)
                self.error_list.append(error)
                self.alpha -= self.eta*(np.matmul(-X.T, y - np.matmul(X,self.alpha)) + self.C*self.alpha)
            if show_result == True:
                sns.set()
                plt.figure()
                plt.plot(self.error_list,c = 'b',label = '$loss quad curve$')
                plt.title('$Gradientquad Descent quad error$')
                plt.xlabel('$epoch$')
                plt.ylabel('$error$')
                plt.legend()
                plt.show()
            return self.alpha,self.error_list
    
        def AdamOptimizer(self,dataset,max_epoch,eta = 0.001,show_result = False):
            """
            parameters:
                input:
                    epsilon     type:float      the stop condition      Default:1e-6 
                    eta    type:float   learning rate   Default:0.001
                    dataset    type:DataFrame    the last columns is label
                    show_result    type:bool    Default:False
                output:
                    if show_result == True:
                        output the figure of the process of training
                return : the parameters of linear model,traing_error
                    type: ndarray
            """
            beta1 = 0.9
            beta2 = 0.999
            X = np.asarray(dataset.iloc[:,:-1])
            y = np.asarray(dataset.iloc[:,-1]).reshape(-1,1)
            self.error_list = []
            # sample the initial params from the unifrom distribution 
            self.alpha = np.random.uniform(0,1,size = [dataset.shape[1] - 1,1])
            # print("##########################",self.alpha.shape)
            error = self.score(dataset)
            m = np.zeros(shape = [X.shape[1],1])
            v = np.zeros(shape = [X.shape[1],1])
            t = 1
            while t < max_epoch:
                error = self.score(dataset)
                self.error_list.append(error)
                grad = np.matmul(-X.T, y - np.matmul(X,self.alpha)) + self.C*self.alpha
                m = beta1*m + (1 - beta1)*grad
                v = beta2*v + (1 - beta2)*np.power(grad,2)
                m_hat = m/(1 - beta1**t)
                v_hat = v/(1 - beta2**t)
                self.alpha -= eta/(np.sqrt(v_hat)+1e-8)*m_hat
                t += 1
            
            if show_result == True:
                sns.set()
                plt.figure(dpi = 110)
                plt.plot(self.error_list,c = 'b',label = '$loss quad curve$')
                plt.title('$Adamquad Optimizer quad error$')
                plt.xlabel('$epoch$')
                plt.ylabel('$error$')
                plt.legend()
                plt.show()
            return self.alpha,self.error_list
    
        def trace(self,dataset,max_epoch,C_list):
            """
            parameters:
                dataset     dtype:dataframe     the last columns is label
                C_list      The list of C
                output      The figure of ridge trace
            """
            #保存原始对象的正则化参数
            init_C = self.C
            beta_list = []
            for C in C_list:
                self.C = C
                beta,_ = self.batch_gradient_descent(dataset,max_epoch)
                beta_list.append(beta)
            beta_list = np.asarray(beta_list)
            plt.figure()
            for i in range(beta_list.shape[1]):
                plt.plot(C_list,beta_list[:,i])
            plt.xlabel('alpha')
            plt.ylabel('weights')
            plt.title('Ridge coefficients as a function of the regularization')
            plt.axis('tight')
            plt.show()
            #恢复初始化参数
            self.C = init_C
    
    #Child Class Lasso
    class Lasso(LinearRegression):
        '''
        Lasso Object(child class)
            attribute:
                c:噪声的影响强度,取值[0,1] 默认为0
                shape:生成数据集的形状,默认为None     dtype = [sample_num,feature_num]
                mean:x所服从正态分布的均值向量       默认为None
                cov:x所服从正态分布的协方差矩阵      默认为None
                noise_miu：噪声所服从的正态分布的均值
                noise_sigma:噪声所服从的正态分布的方差
                noise默认服从标准正态分布
                C:正则化项系数 默认为1
            function:
                generate_data(): 按照实例化对象的初始属性生成数据集
                regression_cond():计算Ridege中 X 的条件数
                regression_correlation():可视化矩阵X的相关系数矩阵
                score():通过训练好的model计算测试误差
                predict():利用已训练好的model预测测试样本
                square_R():goodness of fit
                modified_square_R():modified_square_R
                iterate_ridge:Solve lasso by iterative method of solving ridge multiple times
                coordinate_descent():sovle lasso by coordinate descent
        '''
        def __init__(self,c = 0,shape = None,mean = None,cov = None,noise_miu = 0,noise_sigma = 1,C = 1):
            super(Lasso,self).__init__(c,shape,mean,cov,noise_miu,noise_sigma)
            self.C = C
        def iterate_ridge(self,dataset,max_epoch,tol = 1e-4):
            '''
            iterate_ridge:ce
                Solve lasso by iterative method of solving ridge multiple times
                Iterated Ridge Regression using the approximation
                                |w| =~ norm(w,2)/norm(w,1)
            parameters:
                input:
                    1.dataset      type:DataFrame  
                        notation:the last columns is label and the penultimate column is a vector of all ones,but if your label y has been standardized,you can drop the penultimate column.
                    2.max_epoch     type:int
                        notation:The max number of iterate.
                    3.tol       type:float
                        notation:Threshold with coefficient set to 0. The reason why this parameter is necessary is because RidgeRegression does not have the sparsity.So if we do not set this parameter.The coefficients will close to zero but not equal to it.
                return: lasso estimation result
            '''
            X = dataset.iloc[:,:-1]
            y = dataset.iloc[:,-1]
            # sample the initial params from the unifrom distribution 
            self.alpha = np.random.uniform(0,1,size = dataset.shape[1] - 1)
            for epoch in range(max_epoch):
                temp_matrix = self.C*np.linalg.pinv(np.diag(np.abs(self.alpha)),hermitian = True)
                self.alpha = np.matmul(np.linalg.inv(np.dot(X.T,X) + temp_matrix),np.matmul(X.T,y))
            # 置0时，要避开偏置项
            for i in range(self.alpha.shape[0] - 1):
                if self.alpha[i] < tol:
                    self.alpha[i] = 0
            return self.alpha
    
        def coordinate_descent(self,dataset,eps = 1e-6):
            """
            coordinate_descent:
                Solve lasso by corrdinate descent.
                parameters:
                input:
                    1.dataset      type:DataFrame  
                        notation:the last columns is label and the penultimate column is a vector of all ones,but if your label y has been standardized,you can drop the penultimate column.
                    2.eps           type:float      defalut:1e-6
                        notation:if ||alpha_t - alpha_{t+1} || < eps,we think it is convergence
                return
                    lasso estimation result
            """
            X = dataset.iloc[:,:-1]
            y = dataset.iloc[:,-1]
            n = X.shape[0]
            #初始化alpha
            self.alpha = np.random.uniform(0,1,size = X.shape[1])
            #声明变量，用于浅拷贝alpha，以便于对比前一次的参数与后一次的参数，使用ones的原因是保证第一次可以进入循环
            alpha_copy = -np.ones(X.shape[1])
            # 前一次迭代的参数与后一次迭代的参数差的范数小于eps,则终止迭代
            while(np.linalg.norm(alpha_copy - self.alpha) > eps):
                #浅拷贝alpha
                alpha_copy = self.alpha.copy()
                for l in range(X.shape[1]):
                    #取出第l列
                    X_l = X.iloc[:,l]
                    #将去除第l列的X记为X_drop
                    X_drop = X.drop(X.columns[l],inplace = False,axis = 1)
                    #将去除第l列的alpha记为alpha_drop
                    alpha_drop = np.delete(self.alpha,obj = l)
                    #构建去除第l列后的残差向量
                    r = y - np.dot(X_drop,alpha_drop).reshape(y.shape)
                    #构建与C比较的量
                    compare = (2/n)*np.dot(X_l.T,r)
                    if compare > self.C:
                        self.alpha[l] = n*(compare - self.C)/(2*(X_l**2).sum())
                    elif compare < -self.C:
                        self.alpha[l] = n*(compare + self.C)/(2*(X_l**2).sum())
                    else:
                        self.alpha[l] = 0
            return self.alpha

**测试代码(lasso test)**

from Multiple_linear_regression import Lasso
    import numpy as np
    import pandas as pd
    
    shape = [200,8]
    mean = [1,2,3,4,5,6,7,8]
    cov = np.diag(np.ones((8)))
    
    # beta[-1] 为 bias
    beta = [1,0,3,0,2,0,2,0,0]
    # print(mean,'n',cov)
    # 实例化
    # lasso = Lasso(c = 0.3,shape = shape,mean = mean,cov = cov,noise_miu = 0,noise_sigma = 1,C = 1)
    # 生成数据集
    # dataset = lasso.generate_data(beta,random_state = 42)
    
    # print(dataset.head())
    
    # alpha = lasso.iterate_ridge(dataset,max_epoch = 500,tol = 5e-2)
    
    # X = np.asarray(dataset.iloc[:,:-1])
    # y = np.asarray(dataset.iloc[:,-1])
    
    # print("the estimate of lasso:",alpha)
    # print("the error of lasso:",lasso.score(dataset))
    
    
    # from Multiple_linear_regression import LinearRegression
    
    # lr = LinearRegression()
    # beta,_ = lr.fit(dataset)
    # print("the estimate of OLS:",beta.reshape(-1,))
    # print("the error of linear regression:",((y - np.dot(X,list(beta)).reshape(y.shape))**2).sum())
    
    
    
    # lasso = Lasso(c = 0.3,shape = shape,mean = mean,cov = cov,noise_miu = 0,noise_sigma = 1,C = 1)
    # dataset = lasso.generate_data(beta,random_state = 42)
    # print(dataset.head())
    # alpha = lasso.coordinate_descent(dataset,eps = 1e-6)
    
    # print(alpha)
    # print(lasso.score(dataset))
    
    # import matplotlib.pyplot as plt
    # plt.figure()
    # plt.plot(dataset.iloc[:,-1],label = '$y_{true}$')
    # plt.plot(lasso.predict(dataset.iloc[:,:-1]),label = "$y_{predict}$")
    # plt.legend()
    # plt.show()
    
    
    from sklearn.datasets import load_boston
    from Multiple_linear_regression import Lasso
    import pandas as pd
    import numpy as np
    boston = load_boston()
    X = boston.data
    col_names = boston.feature_names
    X = pd.DataFrame(X,columns = col_names)
    # scaler
    X_norm = (X - X.mean())/(np.sqrt((X.shape[0]-1)/X.shape[0])*X.std())
    X_norm.head()
    y = boston.target
    y = pd.DataFrame({'price':y})
    one_col = pd.DataFrame(np.ones(shape = (X_norm.shape[0],1)))
    one_col.columns = ['one']
    X = pd.concat([X_norm,one_col],axis = 1)
    df = pd.concat([X,y],axis = 1)
    
    print(df.head())
    
    lasso = Lasso(C = 1)
    alpha = lasso.coordinate_descent(df)
    print(alpha)
    print(lasso.score(df))
    
    import matplotlib.pyplot as plt
    plt.figure()
    plt.plot(df.iloc[:,-1],label = '$y_{true}$')
    plt.plot(lasso.predict(df.iloc[:,:-1]),label = "$y_{predict}$")
    plt.legend()
    plt.show()

**测试代码(Multiple\_test)**

import numpy as np
    import pandas as pd
    import matplotlib.pyplot as plt
    from Multiple_linear_regression import LinearRegression
    from Multiple_linear_regression import RidgeRegression
    from Multiple_linear_regression import Lasso
    
    
    # README
    # help(LinearRegression)
    
    shape = [100,5]
    mean = [1,2,3,4,5]
    #正态分布不相关与独立等价，设置cov = I 相当于说明生成的正态分布随机变量彼此独立,理论上来说，这样生成的特征间彼此间共线性很弱
    cov = np.diag(np.ones(5))
    # print(cov.shape[0])
    # beta[-1] 为 bias
    beta = [1,2,3,4,5,1]
    # print(mean,'n',cov)
    # 实例化
    lr = LinearRegression(c = 0.1,shape = shape,mean = mean,cov = cov,noise_miu = 0,noise_sigma = 1)
    # 生成数据集
    dataset = lr.generate_data(beta)
    # print(dataset)
    # 计算线性回归的条件数，值为60左右,说明特征间不存在明显的共线性
    print("***cond:",lr.cond(dataset.iloc[:,:-1]))
    # 生成相关系数矩阵反映特征间的共线情况，得到与上相同的结论
    print(lr.correlation(dataset.iloc[:,:-1]))
    # 训练
    alpha,training_error = lr.fit(dataset)
    # 打印估计出的参数和训练误差
    print('LinearRegression Estimated parameter value：',alpha)
    print('training_error:',training_error)
    
    # 利用sklearn划分训练集和测试集测试score
    from sklearn.model_selection import train_test_split
    Xtrain,Xtest,ytrain,ytest = train_test_split(dataset.iloc[:,:-1],
                                                    dataset.iloc[:,-1],
                                                    test_size = 0.2,
                                                    random_state = 42
                                                    )
    training_data = pd.DataFrame(np.hstack([Xtrain,np.asarray(ytrain).reshape(-1,1)]))
    testing_data = pd.DataFrame(np.hstack([Xtest,np.asarray(ytest).reshape(-1,1)]))
    
    alpha,training_error = lr.fit(training_data)
    print('LinearRegression Estimated parameter value：',alpha)
    print('LinearRegression training_error:',training_error)
    _,testing_error = lr.score(testing_data)
    print('LinearRegression testing_error:',testing_error)
    
    y_predict = lr.predict(Xtest)
    print(y_predict)
    
    
    #RidgeRegression test
    # 实例化对象
    ridge = RidgeRegression(c = 0.1,shape = shape,mean = mean,cov = cov,noise_miu = 0,noise_sigma = 1)
    # 生成数据
    dataset = ridge.generate_data(beta)
    #划分训练集测试集
    from sklearn.model_selection import train_test_split
    Xtrain,Xtest,ytrain,ytest = train_test_split(dataset.iloc[:,:-1],
                                                    dataset.iloc[:,-1],
                                                    test_size = 0.2,
                                                    random_state = 42
                                                    )
    training_data = pd.DataFrame(np.hstack([Xtrain,np.asarray(ytrain).reshape(-1,1)]))
    testing_data = pd.DataFrame(np.hstack([Xtest,np.asarray(ytest).reshape(-1,1)]))
    
    #test
    alpha,training_error = ridge.fit(training_data)
    print('RidgeRegression Estimated parameter value：',alpha)
    print('RidgeRegression training_error:',training_error)
    testing_error,_ = ridge.score(testing_data)
    print('RidgeRegression testing_error:',testing_error)
    
    y_predict = ridge.predict(Xtest)
    print(y_predict)
    
    # 结论，数据非病态下,RidgeRegression与LinearRegression表现差异不大
    
    # 线性回归批量梯度下降测试
    # 实例化
    lr = LinearRegression(c = 0.1,shape = shape,mean = mean,cov = cov,noise_miu = 0,noise_sigma = 1)
    # 生成数据集
    dataset = lr.generate_data(beta)
    alpha,training_error_list = lr.batch_gradient_descent(dataset,max_epoch = 70,eta = 0.00001,show_result = True)
    print('LinearRegression Estimated parameter value：',alpha)
    print('LinearRegression min training_error_list:',min(training_error_list))
    
    #Ridge批量梯度下降测试
    ridge = RidgeRegression(c = 0.1,shape = shape,mean = mean,cov = cov,noise_miu = 0,noise_sigma = 1,C = 1)
    # 生成数据集
    dataset = ridge.generate_data(beta)
    alpha,training_error_list = ridge.batch_gradient_descent(dataset,max_epoch = 70,eta = 0.00001,show_result = True)
    print('RidgeRegression Estimated parameter value：',alpha)
    print('RidgeRegression min training_error_list:',min(training_error_list))
    
    
    # AdamOptimizer
    # 非凸简单问题AdamOptimizer的下降速度比batch gradient descent慢
    shape = [100,5]
    mean = [1,2,3,4,5]
    cov = np.diag(np.ones(5))
    beta = [1,2,3,4,5,1]
    lr = LinearRegression(c = 0.1,shape = shape,mean = mean,cov = cov,noise_miu = 0,noise_sigma = 1)
    dataset = lr.generate_data(beta)
    alpha,training_error_list = lr.AdamOptimizer(dataset,eta = 0.1,max_epoch = 1000,show_result = True)
    print(alpha)
    print(min(lr.error_list))
    # 查看valley处的loss变化情况
    import seaborn as sns
    sns.set()
    plt.figure(dpi = 110)
    plt.plot([x for x in range(200,998)],training_error_list[200:-1],c = 'b',label = '$loss quad curve$')
    plt.title('$Adam quad Optimizer quad error$')
    plt.xlabel('$epoch$')
    plt.ylabel('$error$')
    plt.legend()
    plt.show()
    
    
    #Ridge批量梯度下降测试
    ridge = RidgeRegression(c = 0.1,shape = shape,mean = mean,cov = cov,noise_miu = 0,noise_sigma = 1,C = 1)
    # 生成数据集
    dataset = ridge.generate_data(beta)
    alpha,training_error_list = ridge.AdamOptimizer(dataset,max_epoch = 1000,eta = 0.1,show_result = True)
    print('Adaptive moment estimated(RidgeRegression) parameter value：',alpha)
    print('Adaptive moment estimated(RidgeRegression) min training_error_list:',min(training_error_list))
    
    # 查看valley处的loss变化情况
    import seaborn as sns
    sns.set()
    plt.figure(dpi = 110)
    plt.plot([x for x in range(200,998)],training_error_list[200:-1],c = 'b',label = '$loss quad curve$')
    plt.title('$Adam quad Optimizer quad error$')
    plt.xlabel('$epoch$')
    plt.ylabel('$error$')
    plt.legend()
    plt.show()
    
    from sklearn.linear_model import Ridge

[a6fca30027ebe230129bc2ed1109c7be.png]: /images/20221119/c81d1cd34b674787b2e243b245548a90.png