天天看點

周志華《機器學習》課後習題(第七章):貝葉斯分類

7.1 試使用極大似然法估算回瓜資料集 3.0 中前 3 個屬性的類條件機率.

答:

以第一個屬性色澤為例,其值計數如下:

  • 色澤 烏黑 淺白 青綠
  • 好瓜
  • 否 2 4 3
  • 是 4 1 3
周志華《機器學習》課後習題(第七章):貝葉斯分類

其他兩個屬性同理。

7.2* 試證明:條件獨立性假設不成立時,樸素貝葉斯分類器仍有可能産生最優貝葉斯分類器.

周志華《機器學習》課後習題(第七章):貝葉斯分類
周志華《機器學習》課後習題(第七章):貝葉斯分類

隻有在兩者決策邊界之間(淺黃色區域),其分類情況是不同的,在其他區域,樸素貝葉斯分類結果和最優貝葉斯的分類結果是相同的,是以即便屬性之間獨立性假設不成立,樸素貝葉斯在某些條件(本例中就是屬性機率分布在兩者相交區域之外)下任然是最優貝葉斯分類器。

參考:

《On the Optimality of the Simple Bayesian Classifier under Zero-One Loss》

ps.這個例子就是來自該論文,隻做了一點翻譯工作。論文中給出了更全面的理論證明,和樸素貝葉斯産生最優貝葉斯分類的充分必要條件。本打算看完把理論證明也嘗試複述一遍,但能力有限,一方面沒有了解很透徹,另一方面證明過程有點長感覺表達能力有點不大夠用...

7.3 試程式設計實作拉普拉斯修正的樸素貝葉斯分類器,并以西瓜資料集 3.0 為訓練集,對 p.151 "測1" 樣本進行判别.

han1057578619/MachineLearning_Zhouzhihua_ProblemSets

import numpy as np
import pandas as pd
from sklearn.utils.multiclass import type_of_target
from collections import namedtuple


def train_nb(X, y):
    m, n = X.shape
    p1 = (len(y[y == '是']) + 1) / (m + 2)  # 拉普拉斯平滑

    p1_list = []  # 用于儲存正例下各屬性的條件機率
    p0_list = []

    X1 = X[y == '是']
    X0 = X[y == '否']

    m1, _ = X1.shape
    m0, _ = X0.shape

    for i in range(n):
        xi = X.iloc[:, i]
        p_xi = namedtuple(X.columns[i], ['is_continuous', 'conditional_pro'])  # 用于儲存每個變量的情況

        is_continuous = type_of_target(xi) == 'continuous'
        xi1 = X1.iloc[:, i]
        xi0 = X0.iloc[:, i]
        if is_continuous:  # 連續值時,conditional_pro 儲存的就是 [mean, var] 即均值和方差
            xi1_mean = np.mean(xi1)
            xi1_var = np.var(xi1)
            xi0_mean = np.mean(xi0)
            xi0_var = np.var(xi0)

            p1_list.append(p_xi(is_continuous, [xi1_mean, xi1_var]))
            p0_list.append(p_xi(is_continuous, [xi0_mean, xi0_var]))
        else:  # 離散值時直接計算各類别的條件機率
            unique_value = xi.unique()  # 取值情況
            nvalue = len(unique_value)  # 取值個數

            xi1_value_count = pd.value_counts(xi1)[unique_value].fillna(0) + 1  # 計算正樣本中,該屬性每個取值的數量,并且加1,即拉普拉斯平滑
            xi0_value_count = pd.value_counts(xi0)[unique_value].fillna(0) + 1

            p1_list.append(p_xi(is_continuous, np.log(xi1_value_count / (m1 + nvalue))))
            p0_list.append(p_xi(is_continuous, np.log(xi0_value_count / (m0 + nvalue))))

    return p1, p1_list, p0_list


def predict_nb(x, p1, p1_list, p0_list):
    n = len(x)

    x_p1 = np.log(p1)
    x_p0 = np.log(1 - p1)
    for i in range(n):
        p1_xi = p1_list[i]
        p0_xi = p0_list[i]

        if p1_xi.is_continuous:
            mean1, var1 = p1_xi.conditional_pro
            mean0, var0 = p0_xi.conditional_pro
            x_p1 += np.log(1 / (np.sqrt(2 * np.pi) * var1) * np.exp(- (x[i] - mean1) ** 2 / (2 * var1 ** 2)))
            x_p0 += np.log(1 / (np.sqrt(2 * np.pi) * var0) * np.exp(- (x[i] - mean0) ** 2 / (2 * var0 ** 2)))
        else:
            x_p1 += p1_xi.conditional_pro[x[i]]
            x_p0 += p0_xi.conditional_pro[x[i]]

    if x_p1 > x_p0:
        return '是'
    else:
        return '否'


if __name__ == '__main__':
    data_path = r'C:\Users\hanmi\Documents\xiguabook\watermelon3_0_Ch.csv'
    data = pd.read_csv(data_path, index_col=0)

    X = data.iloc[:, :-1]
    y = data.iloc[:, -1]
    p1, p1_list, p0_list = train_nb(X, y)

    x_test = X.iloc[0, :]   # 書中測1 其實就是第一個資料

    print(predict_nb(x_test, p1, p1_list, p0_list))      

這裡代碼很簡單。不怎麼規範。

7.4 實踐中使用式 (7.15)決定分類類别時,若資料的維數非常高,則機率連乘

周志華《機器學習》課後習題(第七章):貝葉斯分類

的結果通常會非常接近于 0 從試述防止下溢的可能方案.而導緻下溢.

這在p153中已經給出答案了。即取對數将連乘轉化為“連加”防止下溢。

周志華《機器學習》課後習題(第七章):貝葉斯分類

7.5 試證明:二分類任務中兩類資料滿足高斯分布且方差相同時,線性判别分析産生貝葉斯最優分類器.

首先看一下貝葉斯最優分類器:在書中p148中解釋了對于最小化分類錯誤率的貝葉斯最優分類器可表示為:

周志華《機器學習》課後習題(第七章):貝葉斯分類
周志華《機器學習》課後習題(第七章):貝葉斯分類
周志華《機器學習》課後習題(第七章):貝葉斯分類

7.6 試程式設計實作 AODE 分類器,并以西瓜資料集 3.0 為訓練集,對 p.151的"測1" 樣本進行判别.

代碼在:

'''
目前僅拿西瓜資料集測試過,運作正常,其他資料未測試
'''
import numpy as np
import pandas as pd
from sklearn.utils.multiclass import type_of_target

class AODE(object):

    def __init__(self, m):
        self.m_hat = m
        self.m = None
        self.n = None
        self.unique_y = None
        self.n_class = None
        self.is_continuous = None
        self.unique_values = None
        self.total_p = None

    def predict(self, X):
        X = np.array(X)
        if self.total_p == None:
            raise Exception('you have to fit first before predict.')

        result = pd.DataFrame(np.zeros((X.shape[0], self.unique_y.shape[0])), columns=self.unique_y)

        for i in self.total_p.keys():
            result += self._spode_predict(X, self.total_p[i], i)

        return self.unique_y[np.argmax(result.values, axis=1)]

    def fit(self, X, y):
        X = np.array(X)
        self.m, self.n = X.shape
        self.unique_y = np.unique(y)
        self.n_class = self.unique_y.size

        # 這裡轉為list, 是因為貌似type_of_target 有bug, 隻有在pd.Series類型的時候才能解析為'continuous',
        # 在這裡轉為array類型後解析為 'unknown'了
        is_continuous = pd.DataFrame(X).apply(lambda x: (type_of_target(x.tolist()) == 'continuous'))
        self.is_continuous = is_continuous

        unique_values = {}  # 離散型字段的各取值
        for i in is_continuous[~is_continuous].index:
            unique_values[i] = np.unique(X[:, i])

        self.unique_values = unique_values

        # 擷取可以作為父節點的屬性索引,這裡在論文中取值為30; 但在西瓜書中由于樣本很少, 所有直接取0就好
        parent_attribute_index = self._get_parent_attribute(X)

        total_p = {}
        for i in parent_attribute_index:
            p = self._spode_fit(X, y, i)
            total_p[i] = p

        self.total_p = total_p

        return self

    def _spode_fit(self, X, y, xi_index):
        p = pd.DataFrame(columns=self.unique_y, index=self.unique_values[xi_index])  # 儲存各取值下的條件機率
        nunique_xi = self.unique_values[xi_index].size  # 目前屬性的取值數量

        pc_xi_denominator = self.m + self.n_class * nunique_xi  # 計算 p(c, xi) 的分母 |D| + N * N_i

        for c in self.unique_y:
            for xi in self.unique_values[xi_index]:
                p_list = []  # 儲存y取值為c, Xi取值為xi下各個條件機率p(xj|c, xi)和先驗機率p(c, xi)

                c_xi = (X[:, xi_index] == xi) & (y == c)
                X_c_xi = X[c_xi, :]  # y 取值 為c, Xi 取值為xi 的所有資料

                pc_xi = (X_c_xi.shape[0] + 1) / pc_xi_denominator  # p(c, xi)

                # 實際上這裡在j = i時, 個人了解應該是可以跳過不計算的,因為p(xi|c, xi) = 1, 在計算中都是一樣的但這裡為了友善實作,就沒有跳過了。

                for j in range(self.n):
                    if self.is_continuous[j]:  # 字段為連續值, 假設服從高斯分布, 儲存均值和方差
                        # 這裡因為樣本太少。有時候會出現X_c_xi為空或者隻有一個資料的情況, 如何是離散值,依然可以計算;
                        # 但是連續值的情況下,np.mean會報warning, 隻有一個資料時,方差為0
                        # 所有這時, 均值和方差以類别樣本來替代。
                        if X_c_xi.shape[0] <= 1:
                            p_list.append([np.mean(X[y == c, j]), np.var(X[y == c, j])])
                        else:
                            p_list.append([np.mean(X_c_xi[:, j]), np.var(X_c_xi[:, j])])

                    else:
                        # 計算 p(xj|c, xi)
                        condi_proba_of_xj = (pd.value_counts(X_c_xi[:, j])[self.unique_values[j]].fillna(0) + 1) / (
                                X_c_xi.shape[0] + self.unique_values[j].size)
                        p_list.append(np.log(condi_proba_of_xj))
                p_list.append(np.log(pc_xi))  # p(c, xi)在最後的位置

                p.loc[xi, c] = p_list

        return p

    def _spode_predict(self, X, p, xi_index):

        assert X.shape[1] == self.n
        xi = X[:, xi_index]
        result = pd.DataFrame(np.zeros((X.shape[0], p.shape[1])), columns=self.unique_y)  # 儲存每個樣本為不同類别的對數機率值
        for value in p.index:  # 為了可以使用pandas的索引功能, 對于要預測的X值, 每一次循環求同一取值下樣本的條件機率和
            xi_value = xi == value
            X_split = X[xi_value, :]
            for c in p.columns:
                p_list = p.loc[value, c]  # 儲存p(xj|c, xi) 和 p(c, xi)的清單
                for j in range(self.n):  # 周遊所有的條件機率, 将對應的條件機率相加
                    if self.is_continuous[j]:
                        mean_, var_ = p_list[j]
                        result.loc[xi_value, c] += (
                                -np.log(np.sqrt(2 * np.pi) * var_) - (X_split[:, j] - mean_) ** 2 / (2 * var_ ** 2))
                    else:
                        result.loc[xi_value, c] += p_list[j][X_split[:, j]].values

                result.loc[xi_value, c] += p_list[-1]  # 最後加上p(c, xi)

        return result

    def _get_parent_attribute(self, X):
        '''
        基于屬性下各取值的樣本數量,決定哪些屬性可以作為父屬性。
        關于連續值的處理,在《機器學習》書中也沒有提及,AODE的原論文也沒有提及如何處理連續值,
        考慮到若将連續值x_j作為父屬性時,如何求解p(x_i|c, x_j)條件機率會比較麻煩(可以通過貝葉斯公式求解),
        此外AODE本身就是要将取值樣本數量低于m的屬性去除的,從這個角度來說,連續值就不能作為父屬性了。
        是以這裡連續值不作為父屬性
        :param X:
        :return:
        '''

        enough_quantity = pd.DataFrame(X).apply(
            lambda x: (type_of_target(x.tolist()) != 'continuous') & (pd.value_counts(x) > self.m_hat).all())
        return enough_quantity[enough_quantity].index.tolist()


if __name__ == '__main__':
    data_path = r'C:\Users\hanmi\Documents\xiguabook\watermelon3_0_Ch.csv'
    data = pd.read_csv(data_path, index_col=0)

    X = data.iloc[:, :-1]
    y = data.iloc[:, -1]

    aode = AODE(0)
    print(aode.fit(X, y).predict(X.iloc[[0], :]))      

提一下關于連續值處理的問題。這個書中和原論文(好像)都沒有送出,是以按照自己的了解來處理了。考慮到以下,實作過程中不将連續值作為父屬性。

此外AODE本身就是要将取值樣本數量低于一定門檻值(論文中給出的是30)的屬性去除的,從這個角度來說,連續值就不能作為父屬性了,目前其實可以按照區間劃分将連續值離散化。

另外,雖然在樣本這麼小的情況下,看預測結果實際意義不大,但相比于樸素貝葉斯,AODE對于西瓜資料集的拟合更好(錯誤率更低)。

ps.書中給出的式(7.24)有錯誤的,分母的

周志華《機器學習》課後習題(第七章):貝葉斯分類

改正為

周志華《機器學習》課後習題(第七章):貝葉斯分類

,在第十次印刷的時候糾正了,看舊版書的同學要注意了。

7.7 給定 d 個二值屬性的二分類任務,假設對于任何先驗機率項的估算至少需 30 個樣例,則在樸素貝葉斯分類器式 (7.15) 中估算先驗機率項

周志華《機器學習》課後習題(第七章):貝葉斯分類

需 30 x 2 = 60 個樣例.試估計在 AODE 式 (7.23) 中估算先驗機率項

周志華《機器學習》課後習題(第七章):貝葉斯分類

所需的樣例數(分别考慮最好和最壞情形) .

這裡“假設對于任何先驗機率項的估算至少需 30 個樣例”意味着在所有樣本中, 任意

周志華《機器學習》課後習題(第七章):貝葉斯分類

的組合至少出現30次。

周志華《機器學習》課後習題(第七章):貝葉斯分類

這個問題主要基于書中式7.26,就很容易了解了

首先考慮同父結構,根據式7.26,其聯合分布可以表示為:

周志華《機器學習》課後習題(第七章):貝葉斯分類
周志華《機器學習》課後習題(第七章):貝葉斯分類

好久沒更新...罪過,堕落了...前面八題一個月之前就寫好了,一直在看7.7(主要還是懶.)閱讀材料給出的貝葉斯網相關論文(主要是《A Tutorial on Learning With Bayesian Networks》),下面兩題應該還是需要寫代碼實作的,回頭有時間再補把。

7.9 以西瓜資料集 2.0 為訓練集,試基于 BIC 準則建構一個貝葉斯網.

關于貝葉斯網結構的建構,書中p160隻稍微提了一下,不過還是挺好了解的,《A Tutorial on Learning With Bayesian Networks》11節給出了更詳細的描述。比較簡單是方法就是貪心法:

  • 1) 初始化一個網絡結構;
  • 2) 使用E表示目前合法的改變一條邊的操作集合,比如若兩個節點已經有連接配接,那麼合法操作可以删除或者逆轉,如沒有連接配接則可以增加一條邊,目前必須是在保證不會形成回路的情況;
  • 3) 從中選擇一個使得BIC下降最大的一個作為實際操作;
  • 4) 循環2,3直到BIC不再下降。

論文中也給出了其他算法。

有時間再補代碼吧。

繼續閱讀