---
# System prepended metadata

title: Lesson 10 | 簡易數據分析

---

# Lesson 10 | 簡易數據分析

## 第一節：描述性統計(1)

* 這節課我們終於要開始進行數據分析了，我們先簡單介紹幾個描述性統計的函數

    * 請至[這裡](https://linchin.ndmctsgh.edu.tw/data/Example_data.csv)下載本週的範例資料，另外，請同學們先幫我安裝```pandas```以及```numpy```。

        ```python=
        # 載入套件
        import pandas as pd
        import numpy as np
        ```

        ```python=
        data = pd.read_csv("Example_data.csv")
        data
        ```
        
    * 以下幾個```numpy```的函式為常見的描述性統計，分別是：
    
        * 函式「.mean()」可以幫助我們計算算數平均值。
        
            ```python=
            np.mean(data['eGFR'][pd.isna(data['eGFR']) == False])
            ```
            
        * 函數「.std()」可以幫助我們計算標準差。

            ```python=
            np.std(data['eGFR'][pd.isna(data['eGFR']) == False])
            ```
            
        * 函數「.var()」可以幫助我們計算變異數。

            ```python=
            np.var(data['eGFR'][pd.isna(data['eGFR']) == False])
            ```

        * 函數「.median()」可以幫助我們計算中位數。
            
            ```python=
            np.median(data['eGFR'][pd.isna(data['eGFR']) == False])
            ```

        * 函數「quantile()」可以幫助我們計算百分位數。
            
            ```python=
            np.quantile(data['eGFR'][pd.isna(data['eGFR']) == False], [0, 0.25, 0.50, 0.75, 1])
            ```

            ```python=
            # 50百分位數
            np.quantile(data['eGFR'][pd.isna(data['eGFR']) == False], 0.50)
            ```
            
            ```python=
            # 95百分位數
            np.quantile(data['eGFR'][pd.isna(data['eGFR']) == False], 0.95)
            ```            

        * 函數「min()」可以幫助我們找出最小值。
            
            ```python=
            np.min(data['eGFR'][pd.isna(data['eGFR']) == False])
            ```

        * 函數「max()」可以幫助我們找出最大值。

            ```python=
            np.max(data['eGFR'][pd.isna(data['eGFR']) == False])
            ```

        * 函數「table()」可以幫助我們計數各類別的頻率。

            ```python=
            np.unique(data['Disease'], return_counts = True)
            ```
            
        * 稍微加工一下就可以產生各類別的百分比。
            
            ```python=
            np.unique(data['Disease'], return_counts = True)[1]/len(data['Disease'])
            ```
            
## 第一節：描述性統計(2)

* 利用索引，我們將能得到某種條件下的描述性統計，舉例來說我們想獲得Case組的eGFR平均數，可以使用下列語法。

    ```python=
    clean_data = data[(pd.isna(data['eGFR']) == False) & (pd.isna(data['Disease']) == False)]
    np.mean(clean_data['eGFR'][clean_data['Disease'] == 1])
    ```
    
* Control組的亦同

    ```python=
    clean_data = data[(pd.isna(data['eGFR']) == False) & (pd.isna(data['Disease']) == False)]
    np.mean(clean_data['eGFR'][clean_data['Disease'] == 0])
    ```
    
* 如果我們想要把平均數±標準差表示出來，在Python裡的String可以直接使用```+```號進行處理：

    ```python=
    clean_data = data[(pd.isna(data['eGFR']) == False) & (pd.isna(data['Disease']) == False)]
    str(np.mean(clean_data['eGFR'][clean_data['Disease'] == 1])) + ' ± ' + str(np.std(clean_data['eGFR'][clean_data['Disease'] == 1]))
    ```
    
* 我們發現這樣的呈現相當醜，我們可以使用函數「round()」來指定我們想要的小數點位數：

    ```python=
    clean_data = data[(pd.isna(data['eGFR']) == False) & (pd.isna(data['Disease']) == False)]
    m = np.mean(clean_data['eGFR'][clean_data['Disease'] == 1])
    s = np.std(clean_data['eGFR'][clean_data['Disease'] == 1])
    str(round(m, 2)) + ' ± ' + str(round(s, 2))
    ```
    
## 練習1：

* 現在請利用自訂函數功能，寫出兩個函數，規格如下：

    * 函數之input必須為一類別變項及一連續變項。
    
    * 函數中必須設定一個變項，讓使用者能決定顯示小數點的位數。

* 問題1：第一個函數必須輸出在該類別變項為不同類別時，該連續變項之『平均數±標準差』，例如：

    ```python=
    func_1(x = data['Disease'], y = data['eGFR'])
    ```
    
    ```python=
    func_1(x = data['Education'], y = data['SBP'], dig = 1)
    ```
    
    <details>
  <summary>解答</summary>
    
    ```python=
    def func_1(x, y, dig = 3):

        x = x[(pd.isna(x) == False) & (pd.isna(y) == False)]
        y = y[(pd.isna(x) == False) & (pd.isna(y) == False)]

        levels = list(set(x))

        out_list = list()
        for i in levels:
            m = np.mean(y[x == i])
            s = np.std(y[x == i])
            m_s_str = str(round(m, dig)) + ' ± ' + str(round(s, dig))
            out_list.append(m_s_str)

        return out_list
    ```
    
    </details>
    
* 問題2：第二個函數必須輸出在該類別變項為不同類別時，該連續變項之『中位數(25百分位-75百分位)』，例如：

    ```python=
    func_2(x = data['Disease'], y = data['eGFR'])
    ```
    
    ```python=
    func_2(x = data['Education'], y = data['SBP'], dig = 1)
    ```


    <details>
  <summary>解答</summary>
    
    ```python=
    def func_2(x, y, dig = 3):

        x = x[(pd.isna(x) == False) & (pd.isna(y) == False)]
        y = y[(pd.isna(x) == False) & (pd.isna(y) == False)]

        levels = list(set(x))

        out_list = list()
        for i in levels:
            q_list = np.quantile(y[x == i], [0.25, 0.50, 0.75])
            q_str = str(round(q_list[1], dig)) + ' (' + str(round(q_list[0], dig)) + '-' + str(round(q_list[2], dig)) + ')'
            out_list.append(q_str)

        return out_list
    ```
    
    </details>
    
## 第二節：線性迴歸(1)

* 線性迴歸(linear regression)是高中已經教過的基本分析技術，目標是從一群數據點中找出$x_i$與$y_i$之間的線性關係，其方程式如下：

    $$\hat{y_i} = b_{0} + b_{1}x_i$$

* 這裡在高中上課的時候有提到我們要如何解這個$b_0$與$b_1$的係數，我們是透過先定義誤差為$ϵ$：
    
    $$y_i = b_{0} + b_{1}x_i + \epsilon$$

    $$\epsilon = y_i - b_{0} - b_{1}x_i$$

    * 接著再定義找出一組$b0$與$b1$使$\sum{\epsilon}^2$最小化：

    $$\mbox{min}(\sum{\epsilon}^2) = \mbox{min}(\sum\limits_{i=1}^{n}(y_i - b_{0} - b_{1}x_i)^2)$$
    
## 第二節：線性迴歸(2) 

* 若我們要求解線性回歸，我們可以分別對$b_0$與$b_1$進行偏微分得到公式，求解其係數，我們也可以透過第三方套件來完成這件事情。

* 我們接著將使用```scikit-learn```套件，在Python裡[scikit-learn](https://scikit-learn.org/stable/)是最常用來做數據分析以及機器學習的功能強大套件請先安裝：

![L10-1](https://i.imgur.com/rrRtXbB.png)

* 我們可以在[這裡](https://scikit-learn.org/stable/modules/linear_model.html)看到關於線型回歸的相關說明。
    
    ```python=
    from sklearn import linear_model
    reg = linear_model.LinearRegression()
    reg.fit([[0], [1], [2]], [0, 1, 2])
    ```
    
    ```python=
    print(reg.intercept_)
    print(reg.coef_)
    ```
    

## 第三節：線性迴歸(3)

* 假設我們使用的是剛剛的資料，讓我們看看使用「linear_model」的方法：
    
    * 記得要把```numpy array```的shape改成2維

        ```python=
        SBP = data['SBP'][(pd.isna(data['SBP']) == False) & (pd.isna(data['DBP']) == False)].to_numpy()
        DBP = data['DBP'][(pd.isna(data['DBP']) == False) & (pd.isna(data['SBP']) == False)].to_numpy()

        DBP = DBP.reshape((DBP.shape[0], 1))

        reg = linear_model.LinearRegression()
        reg.fit(DBP, SBP)

        # 讓我們看一下預測出來線性回歸的參數。
        print(reg.intercept_, reg.coef_)
        ```

* 這樣的線性回歸方程式可以放入多個變項。

    * 如果我們要用SBP以及DBP同時預測eGFR，可以用下列程式碼：

        ```python=
        na_idx = (pd.isna(data['SBP']) == False) & \
            (pd.isna(data['DBP']) == False) & \
            (pd.isna(data['eGFR']) == False)

        clean_data = data[na_idx]

        x = clean_data[['SBP', 'DBP']].to_numpy()
        y = clean_data['eGFR'].to_numpy()

        reg = linear_model.LinearRegression()
        reg.fit(x, y)
        
        # 也是一樣，讓我們看一下預測出來線性回歸的參數。
        print(reg.intercept_, reg.coef_)
        ```

## 練習2：

* 現在我們想要使用某病人的SBP以及DBP來預測eGFR，所以希望你能寫出一個自訂函數，讓使用者輸入單筆SBP與DBP的值，並為他計算出eGFR的值。

    * 但是這時候我們遇到一個問題，那就是預測公式鐵定是以這樣的型式表現，因此假設有一個新的個案缺少SBP或DBP其中一項，那他將無法使用這個式子：

        $$\hat{eGFR_i} = b_{0} + b_{1}SBP_i + b_{2}DBP_i$$

    * 因此這個函數需要有一個功能，就是當使用者沒有輸入SBP而只輸入DBP時，你必須用這個式子協助他進行預測(依此類推)：

        $$\hat{eGFR_i} = b_{0} + b_{1}DBP_i$$
    
    ```python=
    print(pred_func(SBP = 100))
    print(pred_func(DBP = 100))
    print(pred_func(SBP = 100, DBP = 100))
    ```
    
    <details>
  <summary>解答</summary>
    
    ```python=
    def pred_func(SBP = None, DBP = None):

        na_idx = (pd.isna(data['SBP']) == False) & \
        (pd.isna(data['DBP']) == False) & \
        (pd.isna(data['eGFR']) == False)

        clean_data = data[na_idx]

        if SBP is None and DBP is None:

            return("Please enter at least one variable.")

        elif SBP is None:

            x = clean_data['DBP'].to_numpy()
            y = clean_data['eGFR'].to_numpy()

            x = x.reshape((x.shape[0], 1))

            reg = linear_model.LinearRegression()
            reg.fit(x, y)

            return reg.intercept_ + reg.coef_[0] * DBP

        elif DBP is None:

            x = clean_data['SBP'].to_numpy()
            y = clean_data['eGFR'].to_numpy()

            x = x.reshape((x.shape[0], 1))

            reg = linear_model.LinearRegression()
            reg.fit(x, y)

            return reg.intercept_ + reg.coef_[0] * SBP

        else:

            x = clean_data[['SBP', 'DBP']].to_numpy()
            y = clean_data['eGFR'].to_numpy()

            reg = linear_model.LinearRegression()
            reg.fit(x, y)

            return reg.intercept_ + reg.coef_[0] * SBP + reg.coef_[1] * DBP
    ```
    
    </details>
    
    
* 同學們是否可以想到，我們其實應該以物件導向的設計來撰寫這個功能？

## 第三節：在線性迴歸中使用類別變項(1)

* 你可能有注意到剛剛的資料集中，其中有些變項是被編碼成0、1以及0、1、2，這樣的變項可以放入線性迴歸中嗎?

* 答案是可以的，對於已經編碼成0、1的變項我們可以直接把他放入線性迴歸，但是0、1、2的變項我們就不能直接放入線性迴歸之中。

* 那我們該怎樣處理這種有3個類別以上的類別變項呢?方法很簡單就是要把他改成這樣的編碼格式(他們稱之為one-hot encoding)：

    ```python=
    x = [0, 1, 2, 1, None, 2]
    np.array([0, 1, 0, 1, None, 0,
              0, 0, 1, 0, None, 1]).reshape((6, 2), order='F')
    ```
    
## 第三節：在線性迴歸中使用類別變項(2)

* 至於要怎樣把他轉換成這個格式呢?我們可以用這種方式：

    ```python=
    cato = list(set(data['Education'][pd.isna(data['Education']) == False]))
    print(cato)

    data['education_1'] = 0
    data['education_2'] = 0
    for i in range(len(data['Education'])):
        if data.loc[i, 'Education'] == cato[1]:
            data.loc[i, 'education_1'] = 1
        elif data.loc[i, 'Education'] == cato[2]:
            data.loc[i, 'education_2'] = 1
        else:
            data.loc[i, 'education_1'] = None
            data.loc[i, 'education_2'] = None

    print(data)
    ```
    
* 當然，你也可以在第三方套件裡找到這個功能的函式：

    ```python=
    pd.get_dummies(data['Education'], dummy_na = True)
    ```
    
* 你應該會注意到「get_dummies」與我們想要的不太一樣(NA變新的一欄)，所以我們還要稍微修正一下他：

    ```python=
    def one_hot_coding(cato_df):

        dummies = pd.get_dummies(cato_df, dummy_na = True)

        nan_col = dummies[np.nan]
        dummies[nan_col == 1] = None

        col_name = dummies.columns[pd.isna(dummies.columns) == False]
        col_name = col_name[1:]

        final_dummies = dummies[col_name]

        return final_dummies
    ```
    
    * 稍微看一下結果。

        ```python=
        one_hot_coding(data['Education'])
        ```
    
* 有了這項技術以後，那就可以把剛剛的Income以及Education都納入回歸之內了：

    * 先整理資料：

    ```python=
    data_edu = one_hot_coding(data['Education'])
    data_income = one_hot_coding(data['Income'])

    data_edu = data_edu.rename(columns = {1.0: 'educa_1', 2.0 : 'educa_2'})
    data_income = data_income.rename(columns = {1.0: 'income_1', 2.0 : 'income_2'})

    data1 = pd.concat([data, data_edu, data_income], axis = 1)
    ```

    ```python=
    from sklearn import linear_model

    # 進行線性回歸模型的建立。
    na_idx = (pd.isna(data1['SBP']) == False) & \
        (pd.isna(data1['DBP']) == False) & \
        (pd.isna(data1['Income']) == False) & \
        (pd.isna(data1['Education']) == False) & \
        (pd.isna(data1['eGFR']) == False)

    clean_data1 = data1[na_idx]

    x_var = ['SBP', 'DBP', 'income_1', 'income_2', 'educa_1', 'educa_2']

    x = clean_data1[x_var].to_numpy()
    y = clean_data1['eGFR'].to_numpy()

    reg = linear_model.LinearRegression()
    reg.fit(x, y)

    # 也是一樣，讓我們看一下預測出來線性回歸的參數。
    print(reg.intercept_, reg.coef_)
    ```
    
## 總結
                                           
* 有了關於```pandas```以及```numpy```的熟悉使用，再加上```scikit-learn```我們可以輕易地使用python進行資料分析。

* 由於課程設計的限制，我們並沒有完整的介紹Python裡的統計分析功能，在後面的課程將會著重在機器學習以及深度學習等演算法上。

* 現在你應該有能力獨自完成從蒐集資料至資料分析的整個過程，程式語言將能協助你建立pipeline，讓你能隨時看到最新結果！