作业及代码:https://pan.baidu.com/s/1L-Tbo3flzKplAof3fFdD1w 密码:oin0
本次作业的理论部分:吴恩达机器学习(八)聚类与降维(K-Means,PCA)
编程环境:Jupyter Notebook
本章目录
K-means 聚类
案例1: 给定一个二维数据集,使用kmeans进行聚类。数据集:data/ex7data2.mat
1. K-means 计算过程
import numpy as np
import scipy.io as sio
import matplotlib.pyplot as plt
data1 = sio.loadmat('./data/ex7data2.mat')
X = data1['X']
plt.scatter(X[:,0],X[:,1]) # 绘制原始数据
plt.show()
1. 获取样本点所属类别
计算每个样本点与聚类中心的距离:
c
(
i
)
:
=
j
that minimizes
∥
x
(
i
)
−
μ
j
∥
2
c^{(i)}:=j \quad \text { that minimizes } \quad\left\|x^{(i)}-\mu_{j}\right\|^{2}
c(i):=j that minimizes ∥∥∥x(i)−μj∥∥∥2
def find_centroids(X,centros):
idx = [] #存放样本所属类别索引,即 0,1,2
for i in range(len(X)):
# (2,) (k,2) -> (k,2)
dist =np.linalg.norm((X[i] - centros),axis=1) #(k,) # 每个样本点与聚类中心计算距离(axis=1表示按列计算)
id_i = np.argmin(dist) # 距离排序,取最小距离在数组中的索引号(0,1,2),即将样本点分给对应聚类中心
idx.append(id_i)
return np.array(idx)# 返回值为每个样本点的所属类别索引
init_centros= np.array([[3, 3], [6, 2], [8, 5]]) # 人为给定(随机初始化)3个初始聚类中心
idx = find_centroids(X,centros) # 每个样本点都找到所属类别
#绘制首次聚类结果
plt.scatter(X[:,0],X[:,1],c=idx,cmap='rainbow')#数据集根据索引idx来区分颜色
plt.scatter(init_centros[:,0],init_centros[:,1],
s=200,marker='x',c=[0,1,2],cmap='rainbow')#末尾参数表示符号大小,形状,颜色
2. 计算重心点
def compute_centros(X,idx,k):
centros = [] # 用于存放聚类中心
for i in range(k):
centros_i = np.mean(X[idx == i],axis=0)
centros.append(centros_i)
return np.array(centros)
centros = compute_centros(X,idx,k=3) # 三个聚类中心由初始给定值 变成 当前的平均值
plt.scatter(X[:,0],X[:,1],c=idx,cmap='rainbow')#数据集根据索引idx来区分颜色
plt.scatter(init_centros[:,0],init_centros[:,1],
s=200,marker='x',c=[0,1,2],cmap='rainbow')#初始聚类中心
plt.scatter(centros[:,0],centros[:,1],s=200,marker='x',c='k')# 新的聚类中心
黑色X即为初始聚类的重心点:
3. 运行kmeans,重复执行1和2
def run_kmeans(X,centros,iters):
k = len(centros)
centros_all = []
centros_all.append(centros)# 记录初始聚类中心
centros_i = centros # 初始聚类中心
for i in range(iters):
idx = find_centroids(X,centros_i) # 获取样本点所属类别
centros_i = compute_centros(X,idx,k) # 更新聚类中心
centros_all.append(centros_i) # 记录聚类中心
return idx,np.array(centros_all)
4. 绘制数据集聚类结果&聚类中心的移动轨迹
def plot_data(X,centros_all,idx):
plt.figure()
plt.scatter(X[:,0],X[:,1],c=idx,cmap='rainbow')#数据集根据索引idx来区分颜色
# centros_all是3维数组,三个维度分别是迭代次数,类别索引,特征
plt.plot(centros_all[:,:,0], #横轴
centros_all[:,:,1], #纵轴
'kx--') # k表示黑色,x表示聚类中心,--表示连接线
idx,centros_all = run_kmeans(X,init_centros,iters=10)
plot_data(X,centros_all,idx)
plt.scatter(centros_all[10,:,0],centros_all[10,:,1],s=100,marker='o',c='k')#最终的聚类中心
2. 观察初始聚类点的位置对聚类效果的影响
# 在训练样本 X-中随机选取 k个样本作为初始聚类中心
def init_centros(X,k):
index = np.random.choice(len(X),k)
return X[index]
init_centros(X,k=3)
for i in range(4):
idx,centros_all = run_kmeans(X,init_centros(X,k=3),iters=10)
plot_data(X,centros_all,idx)
plt.scatter(centros_all[10,:,0],centros_all[10,:,1],s=200,marker='o',c='k')#最终的聚类中心
由图可知,不同的初始化位置会得到不同的聚类结果(黑色●),且第一种情况陷入了局部最优。
3. 使用 k-means 压缩图像
-
原始图像:24位真彩色图像,每个像素有RGB三个通道,每个通道有256个强度值(每个像素可以有 2 8 ∗ 2 8 ∗ 2 8 = 2 24 2^8*2^8*2^8=2^{24} 28∗28∗28=224 种色彩)。
-
压缩目标:将每个像素通道的强度值压缩为16种
通俗的讲,原图像每个像素点的单个通道有256种颜色可供选择,而我们要用K-means算法选16种颜色,用于图片压缩。把原始图片的每个像素看作一个数据样本,然后利用K-means算法去找分组最好的16种颜色。
# 原始图像
from skimage import io
image = io.imread('data/bird_small.png')
plt.imshow(image)
# 原始图像数据集
data = sio.loadmat('data/bird_small.mat')
>>> data.keys()
> dict_keys(['__header__', '__version__', '__globals__', 'A'])
A = data['A']
>>> A.shape
> (128, 128, 3) #128×128个像素点,每个像素点有3个通道强度值(范围是0-255)
# 特征归一化
A = A / 255
A = A.reshape(-1,3) #将A转换为二维数组,且列数为3
>>> A.shape # 128*128=16348行(像素点数),3列(通道数)
> (16384, 3)
k = 16
idx,centros_all = run_kmeans(A,init_centros(A,k=16),iters=20)获取原始数据的16个聚类中心
centros = centros_all[-1]#最后的聚类中心
im = np.zeros(A.shape)
for i in range(k):
im[idx==i] = centros[i] # 使同一类点全都相等于聚类中心点
im = im.reshape(128,128,3) # 将聚类后的数据重构为图像格式
plt.imshow(im)
PCA 降维
1. PCA 计算过程
j计算过程请参考 :聚类与降维(K-Means,PCA)# 2.2 节
数据集:data/ex7data1.mat。
要求:将二维数据降到一维。
import numpy as np
import scipy.io as sio
import matplotlib.pyplot as plt
mat = sio.loadmat('./data/ex7data1.mat')
X = mat['X']
>>> X.shape
> (50, 2)
plt.scatter(X[:,0],X[:,1])
plt.show()
1. 对X去均值化,使得特征的数值在可比较的范围之内
X_demean = X - np.mean(X,axis=0)
plt.scatter(X_demean[:,0],X_demean[:,1])
plt.show()
2. 计算协方差矩阵
S
i
g
m
a
=
1
m
X
T
X
Sigma=\frac{1}{m} X^{T} X
Sigma=m1XTX
Sigma = X_demean.T @ X_demean / len(X)
3. 利用奇异值分解(SVD)来分解协方差矩阵,得到特征矩阵
U
n
×
n
U_{n×n}
Un×n
[
U
,
S
,
V
]
=
s
v
d
(
S
i
g
m
a
)
[U, S, V]= svd(Sigma)
[U,S,V]=svd(Sigma)
U
n
×
n
=
[
∣
∣
∣
∣
∣
u
(
1
)
u
(
2
)
u
(
3
)
…
u
(
n
)
∣
∣
∣
∣
∣
]
U_{n \times n}=\left[\begin{array}{ccccc} | & | & | & | & | \\ u^{(1)} & u^{(2)} & u^{(3)} & \ldots & u^{(n)} \\ | & | & | & | & | \end{array}\right]
Un×n=⎣⎡∣u(1)∣∣u(2)∣∣u(3)∣∣…∣∣u(n)∣⎦⎤
U,S,V = np.linalg.svd(Sigma)
>>> U
> array([[-0.76908153, -0.63915068], # 可以发现两个特征向量正交
[-0.63915068, 0.76908153]])
plt.figure(figsize=(7,7))
plt.scatter(X_demean[:,0],X_demean[:,1])
plt.plot([0,U1[0]],[0,U1[1]],c='r') # 红线代表主成分特征向量
plt.plot([0,U[:,1][0]],[0,U[:,1][1]],c='k') # 黑线代表次成分特征向量
plt.show()
4. 实现降维
z
(
i
)
=
U
r
e
d
u
c
e
T
x
(
i
)
=
[
∣
∣
∣
∣
∣
u
(
1
)
u
(
2
)
u
(
3
)
…
u
(
k
)
∣
∣
∣
∣
∣
]
T
x
(
i
)
z^{(i)}=U_{r e d u c e}^{T}\ x^{(i)}=\left[\begin{array}{ccccc} | & | & | & | & | \\ u^{(1)} & u^{(2)} & u^{(3)} & \ldots & u^{(k)} \\ | & | & | & | & | \end{array}\right]^{T}x^{(i)}
z(i)=UreduceT x(i)=⎣⎡∣u(1)∣∣u(2)∣∣u(3)∣∣…∣∣u(k)∣⎦⎤Tx(i)
U_reduce = U[:,0] # 取第一列作为降维后的特征向量 #(2,)
U_reduce = U_reduce.reshape((2,1)) #(2,1)
z = X_demean @ U_reduce
>>> z.shape
> (50,1) # 可以看出维度由(50,2)降至(50,1)
5. 还原数据
x
a
p
p
r
o
x
=
U
r
e
d
u
c
e
⋅
z
≈
x
x_{a p pr o x}=U_{r e d u c e} \cdot z\approx x
xapprox=Ureduce⋅z≈x
X_approx = z @ U_reduce.T + np.mean(X,axis=0) #(50,2) 加上均值是为了还原到原始坐标系
plt.scatter(X[:,0],X[:,1])
plt.scatter(X_approx[:,0],X_approx[:,1])
2. 使用 PCA 对图像降维
数据集:data/ex7faces.mat,要求将特征降至 k = 36 维
import numpy as np
import scipy.io as sio
import matplotlib.pyplot as plt
mat = sio.loadmat('data/ex7faces.mat')
X = mat['X']
>>> X.shape
> (5000, 1024)
可视化前100个人脸原始图像:
def plot_100_images(X):
fig, axs = plt.subplots(ncols=10, nrows=10, figsize=(10,10))
for c in range(10):
for r in range(10):
axs[c,r].imshow(X[10*c + r].reshape(32,32).T,cmap = 'Greys_r') #显示单通道的灰度图
axs[c,r].set_xticks([])
axs[c,r].set_yticks([])
plot_100_images(X)
means = np.mean(X,axis=0)
X_demean = X - means
Sigma = X_demean.T @ X_demean / len(X)
U,S,V = np.linalg.svd(Sigma)
U_reduce = U[:,:36]
z = X_demean @ U_reduce
>>> z.shape
> (5000, 36)
X_approx = z @ U_reduce.T
plot_100_images(X_appox)