python 数据分析 day3 矩阵及利用线性回归处理boston数据集

Python 数据分析day_1:认识Notebook，及NumPy , Pandas 和 Matplotlib基本使用方法
python 数据分析 day2 数组的运算及矩阵
python 数据分析 day3 矩阵及利用线性回归处理boston数据集
python 数据分析 day4 pandas 之 Series
python 数据分析day4 Pandas 之 DataFrame

矩阵

定义矩阵np.matrix();推荐np.array()方法创建矩阵

m1=np.matrix('1 2 3;4 5 6')
"""
matrix([[1, 2, 3],
        [4, 5, 6]])
"""

矩阵的属性

矩阵对应的数组

m1.A
"""
array([[1, 2, 3],
       [4, 5, 6]])
"""

矩阵对应的一维数组

m1.A1 #array([1, 2, 3, 4, 5, 6])
"""
"""

逆矩阵

m1.I
"""
matrix([[-0.94444444,  0.44444444],
        [-0.11111111,  0.11111111],
        [ 0.72222222, -0.22222222]])
"""

矩阵转置:X.T/X.transpose/X.swapaxes(0,1)

m1.T
"""
matrix([[1, 4],
        [2, 5],
        [3, 6]])
"""
m1.transpose()
"""
matrix([[1, 4],
        [2, 5],
        [3, 6]])
"""
m1.swapaxes(1,0)
matrix([[1, 4],
        [2, 5],
        [3, 6]])

矩阵的运算

m1
"""
matrix([[1, 2, 3],
        [4, 5, 6]])
"""
m2=np.array([[1,2],[3,4],[5,6]])
"""
matrix([[1, 2],
        [3, 4],
        [5, 6]])
"""

矩阵乘法:m1*m2/m1@m2/数乘/np.matmul()

m1*m2
"""
matrix([[22, 28],
        [49, 64]])
"""
m1 @ m2
"""
matrix([[22, 28],
        [49, 64]])
"""
m1*3
"""
array([[ 3,  6,  9],
       [12, 15, 18]])
"""
np.matmul(m1,m2)
"""
matrix([[22, 28],
        [49, 64]])
"""

求行列式的值:np.linalg.dat()

array3 = np.array([[1, 2], [3, 4]])

"""
array([[1, 2],
       [3, 4]])
"""
np.linalg.det(array3) #-2.0000000000000004

矩阵对应的行列式的值是0，矩阵是一个奇异矩阵，无法求逆

array5 = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
"""
array([[1, 2, 3],
       [4, 5, 6],
       [7, 8, 9]])
"""
np.linalg.det(array5) #0.0

求矩阵的秩：np.linalg.matrix_rank(X)

array6 = np.array([[1, 1, 1], [2, 2, 2], [3, 3, 3]])
np.linalg.matrix_rank(array6) #1

求解方程组：

方法一：np.linalg.solve(X,y)

$$\{\begin{array}{l}x+2y+z=8\\ 3x+7y+2z=23\\ 2x+2y+z=9\end{array}$$

$$Ax=b$$

$$rank(A)=rank(Ab)$$

A = np.array([[1,2,1],[3,7,2],[2,2,1]])
"""
A
array([[1, 2, 1],
       [3, 7, 2],
       [2, 2, 1]])
"""

np.linalg.matrix_rank(A)
# 3

b = np.array([8,23,9]).reshape(-1,1)
"""
b
array([[ 8],
       [23],
       [ 9]])
"""

合成增广矩阵

np.hstack((A,b))
"""
array([[ 1,  2,  1,  8],
       [ 3,  7,  2, 23],
       [ 2,  2,  1,  9]])
"""

计算增广矩阵的秩

np.linalg.matrix_rank(np.hstack((A,b)))
# 3

求解

np.linalg.solve(A,b)
"""
array([[1.],
       [2.],
       [3.]])
"""

法二：利用逆矩阵直接带公式求解

np.linalg.inv(A) @ b
"""
array([[1.],
       [2.],
       [3.]])
"""

矩阵堆叠：np.vstack((x,y))/np.hstack((x,y))

# 垂直堆叠
x = np.array([[1, 1, 1], [2, 2, 2], [3, 3, 3]])
y = np.array([[4, 4, 4], [5, 5, 5], [6, 6, 6]])
np.vstack((x, y))

"""
array([[1, 1, 1],
       [2, 2, 2],
       [3, 3, 3],
       [4, 4, 4],
       [5, 5, 5],
       [6, 6, 6]])
"""
# 水平堆叠
np.hstack((x, y))
"""
array([[1, 1, 1, 4, 4, 4],
       [2, 2, 2, 5, 5, 5],
       [3, 3, 3, 6, 6, 6]])
"""

“”"

“”"

线性回归

获取并处理波士顿房价数据

导入数据及通过内置函数dir()了解对象的属性和方法

from sklearn.datasets import load_boston

dataset = load_boston()
# 通过Python内置函数dir了解对象的属性和方法
dir(dataset) #['DESCR', 'data', 'feature_names', 'filename', 'target']

通过pd.DataFrame()将数据以表展示

df = pd.DataFrame(data=dataset.data, columns=dataset.feature_names)
df

将房价数据拿出来，也放到表里:df[‘PRICE’]=dataset.target

df['PRICE'] = dataset.target
df

计算协方差：df.cov() ,缺点是不太直观，标准化处理之后更直观

计算相关系数：df.corr()，标准化后的协方差，默认给皮尔逊系数

探索房间数和房价的关系

绘制房间数和房价的对应的散点图

x, y = df['RM'], df['PRICE']
plt.scatter(x, y)

将房间数和对应的房价保存到字典中

my_dict = {key: value for key, value in zip(x, y)}
my_dict

处理字典型数据类型的示例

直接对字典排序，排序对象是键，不是值

prices = {
    'AAPL': 191.88,
    'GOOG': 1186.96,
    'IBM': 149.24,
    'ORCL': 48.44,
    'ACN': 166.89,
    'FB': 208.09,
    'SYMC': 21.29
}
sorted(prices) #['AAPL', 'ACN', 'FB', 'GOOG', 'IBM', 'ORCL', 'SYMC']

要想对字典的值进行排序，可以使用匿名函数
下面是以值进行排序后的，生成的对应的键的数组

sorted(prices, key=lambda x: prices[x]) #['SYMC', 'ORCL', 'IBM', 'ACN', 'AAPL', 'FB', 'GOOG']

kNN算法进行预测（k-Nearest-Neighbors）

找到跟当前数据最接近的k条数据，用它们的情况来预估当前的情况。

# 输入小区平均房间数量，返回预测的房价
def knn(history_data, rm, top_k=5):
    most_near_key = sorted(history_data, key=lambda x: (x - rm) ** 2)[:top_k]
    return np.mean([my_dict[x] for x in most_near_key])

knn(my_dict, 6.5, top_k=3) #24.2

knn(my_dict, 7.2, top_k=10) #35.26

损失函数

定义损失函数

为了判定预测效果（拟合效果）是否理想，引入MSE（均方误差, Mean Squared Error）损失函数。

$$MSE=\frac{1}{N}\sum _{i=1}^N(y_i-\widehat{y_i}{)}^2$$

# 找到一条线对图上的点做最佳拟合 ---> 最佳的标准：误差的平方和最小 ---> 损失函数
def get_loss(y, y_hat):
    return np.mean((y - y_hat) ** 2)

蒙特卡洛模拟进行预测

$$y=ax+b$$

# 用随机的方式产生a和b的值，试图让MSE变得更小
min_loss = np.inf
best_a, best_b = None, None
for i in range(10000):
    a, b = np.random.random(2) * 200 - 100
    y_hat = a * x + b
    curr_loss = get_loss(y, y_hat)
    # 找到更好的拟合
    if curr_loss < min_loss:
        min_loss = curr_loss
        best_a, best_b = a, b
        print(f'第{i}次', best_a, best_b, min_loss)

得到a,b的值后，就可计算预测值了，a,b的值每次进行实验都可能不同

y_hat1 = 9.067 * x - 34.233
y_hat1
y_hat2 = 9.3163 * x - 35.531
y_hat2

绘制图像查看拟合情况

plt.scatter(x, y)
plt.plot(x, y_hat1, linewidth=4, color='red')
plt.plot(x, y_hat2, linewidth=4, color='green')
plt.show()

梯度下降法进行预测

#python 梯度下降算法实现
def partial_a(a,b,x,y):
#     return 2*np.mean((y-a*x-b)*(-x))
    return 2*np.mean((a*x+b-y)*x)

def partial_b(a,b,x,y):
#     return 2*np.mean(-y+a*x+b)
    return 2*np.mean((a*x+b-y))

# a,b=np.random.random(2)*200-100
a,b=20,-35
# min_loss=get_loss()
delta=0.01
for _ in range(1000):
    a=a-partial_a(a,b,x,y)*delta
    b=b-partial_b(a,b,x,y)*delta
curr_loss=get_loss(y,a*x+b)
print(a,b,curr_loss) #9.184417171048281 -35.19418978913763 43.603929193890764

使用最小二乘法的内置函数进行预测：np.linalg.lstsq()

# 最小二乘法
A=np.vstack([x,np.ones(len(x))]).T
# print(x,x.shape,type(x))
print(A)
B=y
results=np.linalg.lstsq(A,B.values,rcond=None)
a,b=results[0]
loss=results[1][0]/506
print(a,b)
print(loss)

"""
[[6.575 1.   ]
 [6.421 1.   ]
 [7.185 1.   ]
 ...
 [6.976 1.   ]
 [6.794 1.   ]
 [6.03  1.   ]]
9.102108981180313 -34.67062077643857
43.600551771169584
"""

正规方程法求系数

#正规方程法求系数
X=np.vstack([np.ones(len(x)),x]).T
# print(X)
def normalEqn(X,y):
    theta=np.linalg.inv(X.T@X)@X.T@y
    return theta

print(normalEqn(X,y)) #[-34.67062078   9.10210898]