import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
# 三維，兩個特征
fig = plt.figure(figsize=(8, 6)) #設置圖標的大小
ax = fig.add_subplot(111, projection='3d') # 111的意思是把畫布分為1行1列，畫在第一個方格內(nèi)。其實就是整個畫布。

# 創(chuàng)建樣本，注意兩個特征不能線性相關，否則無法用最小二乘解參數(shù)
X1 = np.arange(-4, 4, 0.1)
m = len(X1)
X2 = np.random.rand(m)*5
# print(X2)
# print(X1)

# 堆疊全1數(shù)組和X1以及X2形成樣本的矩陣，倒置，用以矩陣乘法
X = np.vstack((np.full(m, 1), X1, X2)).T

# y = 15*X1 + 3 * X2 + theta0
# 自定義樣本輸出
Y = X1 + 3 * X2 + 3*np.random.randn(m)

# 利用標準方程(最小二乘法求解theta)
theta = np.dot(np.dot(np.linalg.inv(np.dot(np.transpose(X), X)),
np.transpose(X)), Y)
print(theta)

# 構造網(wǎng)格 meshgrid函數(shù)可以通過延伸矩陣構造多維坐標中的網(wǎng)格坐標。
M, N = np.meshgrid(X1, X2)

# zip函數(shù)構建一個多元祖[(x1,y1),(x2,y2)...],ravel函數(shù)將一個多維數(shù)組合并成一維數(shù)組
Z = np.array([theta[1] * d + theta[2]*p + theta[0] for d, p in zip(np.ravel(M), np.ravel(N))]).reshape(M.shape)

# 根據(jù)網(wǎng)格和函數(shù)構建圖形 suface是一個完整的面
ax.plot_surface(M, N, Z)
# scatter是散點圖
ax.scatter(X1, X2, Y, c='r')
# 設置坐標軸的名稱
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Z')

plt.show()