正文 

当前位置： 首页 > 技术 > java教程

第二课第一周大作业--构建和评估一个线性风险模型(代码片段)

Tina姐 Tina姐     2022-12-04     148

关键词：

之前教程：
第二课第一周第1节-AI用于医学预后简介
第二课第一周第2节-做医学预后，你需要掌握什么？
第二课第一周第3-4节-什么是预后?
第二课第一周第4-7节 医学预后案例欣赏+作业解析
第二课第一周第8节 风险得分计算+作业解析
第二课第一周第9-11节 评估预后模型+作业解析

第一周课程完美结束，让我们来预测糖尿病视网膜病变的风险
作业地址：gongzhonghao 相关文章上领取

outline

大概内容如下

overview

在本作业中，您将使用逻辑回归构建糖尿病患者视网膜病变的风险评分模型。
在开发模型时，我们将了解以下概念：

  • Data preprocessing
    • Log transformations
    • Standardization
  • Basic Risk Models
    • Logistic Regression
    • C-index
    • Interactions Terms

Diabetic Retinopathy

视网膜病变是一种眼部疾病，会导致称为视网膜的眼部血管发生变化。
这通常会导致视力改变或失明。
众所周知，糖尿病患者患视网膜病变的风险很高。

Logistic Regression

逻辑回归是用于预测二元结果概率的适当分析。在我们的案例中，这将是患有或不患有糖尿病视网膜病变的可能性。

逻辑回归是最常用的二元分类算法之一。它用于找到最佳拟合模型来描述一组特征（也称为输入、独立、预测或解释变量）和二元结果标签（也称为输出、依赖或响应）之间的关系。

逻辑回归具有输出预测始终在 [ 0 , 1 ] [0,1] [0,1] 范围内的特性

1 import packages

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

2 load data

from utils import load_data

# This function creates randomly generated data
# X, y = load_data(6000)

# For stability, load data from files that were generated using the load_data
X = pd.read_csv('X_data.csv',index_col=0)
y_df = pd.read_csv('y_data.csv',index_col=0)
y = y_df['y']

这里加载官方提供的数据，如果你没有下载到，也可以使用随机生成数据

X 和 y 是 Pandas DataFrames，包含 6,000 名糖尿病患者的数据。

3 explore the dataset

The features (X) include the following fields:

  • Age: (years)
  • Systolic_BP: Systolic blood pressure (mmHg)收缩压
  • Diastolic_BP: Diastolic blood pressure (mmHg)舒张压
  • Cholesterol: (mg/DL) 胆固醇

target (y) 是患者是否发生视网膜病变的指标

  • y = 1：患者患有视网膜病变。
  • y = 0：患者没有视网膜病变。

在我们建立模型之前，让我们仔细看看我们的训练数据的分布。为此，我们将使用 75/25 拆分将数据拆分为训练集和测试集。

from sklearn.model_selection import train_test_split

X_train_raw, X_test_raw, y_train, y_test = train_test_split(X, y, train_size=0.75, random_state=0)

for col in X.columns:
    X_train_raw.loc[:, col].hist()
    plt.title(col)
    plt.show()

正如我们所看到的，这些分布通常呈钟形分布，但有轻微的向右偏斜。
许多统计模型假设数据是正态分布的，形成一个对称的高斯钟形（没有偏斜），更像下面的例子。

我们可以通过消除偏斜将数据转换为更接近正态分布。消除偏斜的一种方法是将log函数应用于数据。让我们绘制特征变量的log，看看它是否产生了预期的效果。

for col in X_train_raw.columns:
    np.log(X_train_raw.loc[:, col]).hist()
    plt.title(col)
    plt.show()

这里使用了np.log()

对比年龄的分布，取 log 后我们可以看到数据更加对称

4 Mean-Normalize the data

现在让我们 transform 我们的数据，使分布更接近标准正态分布。
首先，我们将使用 log 转换从分布中消除一些偏斜。然后我们将“标准化”分布，使其均值为零，标准差为 1。

编写一个函数，首先去除数据中的一些偏斜，然后标准化分布，以便对于每个数据点 𝑥 𝑥 x ,
x ‾ = x − m e a n ( x ) s t d ( x ) \\overlinex = \\fracx - mean(x)std(x) x=std(x)x−mean(x)​
请注意，我们 test data 是“看不见的”数据。

def make_standard_normal(df_train, df_test):
    """
    In order to make the data closer to a normal distribution, take log
    transforms to reduce the skew.
    Then standardize the distribution with a mean of zero and standard deviation of 1. 
  
    Args:
      df_train (dataframe): unnormalized training data.
      df_test (dataframe): unnormalized test data.
  
    Returns:
      df_train_normalized (dateframe): normalized training data.
      df_test_normalized (dataframe): normalized test data.
    """
    
    ### START CODE HERE (REPLACE INSTANCES OF 'None' with your code) ###  
    # Remove skew by applying the log function to the train set, and to the test set
    df_train_unskewed = np.log(df_train)
    df_test_unskewed = np.log(df_test)
    
    #calculate the mean and standard deviation of the training set
    mean = df_train_unskewed.mean(axis=0)
    stdev = df_train_unskewed.std(axis=0)
    
    # standardize the training set
    df_train_standardized = (df_train_unskewed - mean) / stdev
    
    # standardize the test set (see instructions and hints above)
    df_test_standardized = (df_test_unskewed - mean) / stdev
    
    ### END CODE HERE ###
    return df_train_standardized, df_test_standardized

通过测试可以看到，经过标准化后，训练数据的均值为1，方差为0，测试数据集上虽然不是0和1，但也很接近。

把我们的X_train, X_test进行标准化

X_train, X_test = make_standard_normal(X_train_raw, X_test_raw)

经过数据转换后的数据分布👇

7 build the model

现在我们准备通过使用我们的数据训练逻辑回归来构建风险模型。

def lr_model(X_train, y_train):
    
    ### START CODE HERE (REPLACE INSTANCES OF 'None' with your code) ###
    # import the LogisticRegression class
    from sklearn.linear_model import LogisticRegression
    
    # create the model object
    model = LogisticRegression()
    
    # fit the model to the training data
    model.fit(X_train,y_train)
    
    ### END CODE HERE ###
    #return the fitted model
    return model

model_X = lr_model(X_train, y_train)

8 评估模型

我们使用c-index评估模型

  • c-index 衡量风险评分的辨别力。
  • 直观地说，较高的 c 指数表明模型的预测与一对患者的实际结果一致。
  • c-index的公式是

KaTeX parse error: Undefined control sequence: \\mbox at position 2: \\̲m̲b̲o̲x̲cindex = \\fra…

  • 允许的配对是一对具有不同结果（y=0和y=1）的患者。
  • concordant配对是允许的配对，其中风险评分较高的患者也有较差的结果。
  • ties 是允许配对中，其中患者具有相同的风险评分。
def cindex(y_true, scores):
    '''

    Input:
    y_true (np.array): a 1-D array of true binary outcomes (values of zero or one)
        0: patient does not get the disease
        1: patient does get the disease
    scores (np.array): a 1-D array of corresponding risk scores output by the model

    Output:
    c_index (float): (concordant pairs + 0.5*ties) / number of permissible pairs
    '''
    n = len(y_true)
    assert len(scores) == n

    concordant = 0
    permissible = 0
    ties = 0
    
    ### START CODE HERE (REPLACE INSTANCES OF 'None' with your code) ###
    # use two nested for loops to go through all unique pairs of patients
    for i in range(n):
        for j in range(i+1, n): #choose the range of j so that j>i
            
            # Check if the pair is permissible (the patient outcomes are different)
            if y_true[i] != y_true[j]:
                # Count the pair if it's permissible
                permissible =permissible + 1

                # For permissible pairs, check if they are concordant or are ties

                # check for ties in the score
                if scores[i] == scores[j]:
                    # count the tie
                    ties = ties + 1
                    # if it's a tie, we don't need to check patient outcomes, continue to the top of the for loop.
                    continue

                # case 1: patient i doesn't get the disease, patient j does
                if y_true[i] == 0 and y_true[j] == 1:
                    # Check if patient i has a lower risk score than patient j
                    if scores[i] < scores[j]:
                        # count the concordant pair
                        concordant = concordant + 1
                    # Otherwise if patient i has a higher risk score, it's not a concordant pair.
                    # Already checked for ties earlier

                # case 2: patient i gets the disease, patient j does not
                if y_true[i] == 1 and y_true[j] == 0:
                    # Check if patient i has a higher risk score than patient j
                    if scores[i] > scores[j]:
                        #count the concordant pair
                        concordant = concordant + 1
                    # Otherwise if patient i has a lower risk score, it's not a concordant pair.
                    # We already checked for ties earlier

    # calculate the c-index using the count of permissible pairs, concordant pairs, and tied pairs.
    c_index = (concordant + (0.5 * ties)) / permissible
    ### END CODE HERE ###
    
    return c_index

我们在上节课中，已经学习了c-index的理论基础，这里直接给代码，很好理解。

9 Evaluate the Model on the Test Set

使用 predict_proba 获得测试数据的概率。

scores = model_X.predict_proba(X_test)[:, 1]
c_index_X_test = cindex(y_test.values, scores)
print(f"c-index on test set is c_index_X_test:.4f")

out: c-index on test set is 0.8182

由于线性回归模型，每个特征都有自己的系数（coefficients），也可以理解为每个特征的权重，我们把系数画出来，看哪个特征对结果的影响更大

可以使用函数model.coef_计算系数

coeffs = pd.DataFrame(data = model_X.coef_, columns = X_train.columns)
coeffs.T.plot.bar(legend=None);


从图上可以看出，年龄对视网膜病变的影响最大。

10 Improve the Model

我们可以通过增加交互项改进模型

  • For example, if we have data
    x = [ x 1 , x 2 ] x = [x_1, x_2] x=[x1​,x2​]
  • We could add the product so that:
    x ^ = [ x 1 , x 2 , x 1 ∗ x 2 ] \\hatx = [x_1, x_2, x_1*x_2] x^=[x1​,x2​,x1​∗x2​]

把所有的特征进行两两组合，使用乘法而不是加法（上一节课已经讨论过，为什么用乘法）

把组合后的特征都加入X_train中，重新训练模型。我们得到以下所有特征的系数，加入交互项后，c-index=0.8281,提升了一点一点性能。

您可能会注意到 Age、Systolic_BP 和 Cholesterol 具有正系数。这意味着这三个特征的值越高，疾病的预测概率就越高。您还可能注意到年龄 x 胆固醇的交互作用具有负系数。这意味着年龄 x 胆固醇乘积的较高值会降低疾病的发病概率。

这部分代码略，详情自行下载。

文章持续更新，可以关注微信公众号【医学图像人工智能实战营】获取最新动态，一个关注于医学图像处理领域前沿科技的公众号。坚持已实践为主，手把手带你做项目，打比赛，写论文。凡原创文章皆提供理论讲解，实验代码，实验数据。只有实践才能成长的更快，关注我们，一起学习进步~

我是Tina, 我们下篇博客见~

白天工作晚上写文，呕心沥血

觉得写的不错的话最后，求点赞，评论，收藏。或者一键三连

第二课第一周大作业--构建和评估一个线性风险模型(代码片段)

之前教程：第二课第一周第1节-AI用于医学预后简介第二课第一周第2节-做医学预后，你需要掌握什么？第二课第一周第3-4节-什么是预后?第二课第一周第4-7节医学预后案例欣赏+作业解析第二课第一周第8节风险得分... 查看详情

第二课第一周第9-11节评估预后模型+作业解析

将介绍预测模型的评估。评估预后模型背后的基本思想是查看它在成对患者身上的表现如何。假如这两个苹果为代表的病人。左手的苹果看起来已经不新鲜了，两天后就会过期，但我右手的苹果在接下来的两天里都是新... 查看详情

第二课第一周第9-11节评估预后模型+作业解析

将介绍预测模型的评估。评估预后模型背后的基本思想是查看它在成对患者身上的表现如何。假如这两个苹果为代表的病人。左手的苹果看起来已经不新鲜了，两天后就会过期，但我右手的苹果在接下来的两天里都是新... 查看详情

第二课第一周第9-11节评估预后模型+作业解析

将介绍预测模型的评估。评估预后模型背后的基本思想是查看它在成对患者身上的表现如何。假如这两个苹果为代表的病人。左手的苹果看起来已经不新鲜了，两天后就会过期，但我右手的苹果在接下来的两天里都是新... 查看详情

第二课第一周第8节风险得分计算+作业解析(代码片段)

在本课中，我们将把前面示例中使用风险方程式得到的分数进行形式化地描述，我们还将研究交互项以及它们为何有用。到目前为止，我们看到的预测模型是通过所有特征乘以与该特征相关的系数。例如，年龄特... 查看详情

第二课第一周第8节风险得分计算+作业解析(代码片段)

在本课中，我们将把前面示例中使用风险方程式得到的分数进行形式化地描述，我们还将研究交互项以及它们为何有用。到目前为止，我们看到的预测模型是通过所有特征乘以与该特征相关的系数。例如，年龄特... 查看详情

第二课第一周第4-7节医学预后案例欣赏+作业解析(代码片段)

第二课第一周第4-7节医学预后案例欣赏+作业解析视频地址：B>Tina-姐预后的案例在这节课中，我们将看看预后临床例子。看看预后任务的输入和输出是什么样子的，以及我们如何权衡输入之间的权重。我们可以把... 查看详情

第二课第一周第4-6节医学预后案例欣赏+作业解析(代码片段)

第二课第一周第4-7节医学预后案例欣赏+作业解析视频地址：B>Tina-姐预后的案例在这节课中，我们将看看预后临床例子。看看预后任务的输入和输出是什么样子的，以及我们如何权衡输入之间的权重。我们可以把... 查看详情

第二课第一周第8节风险得分计算+作业解析(代码片段)

在本课中，我们将把前面示例中使用风险方程式得到的分数进行形式化地描述，我们还将研究交互项以及它们为何有用。到目前为止，我们看到的预测模型是通过所有特征乘以与该特征相关的系数。例如，年龄特... 查看详情

第一课第一周大作业-胸部14种疾病分类-代码详解(代码片段)

深度学习胸部X射线诊断本次作业文件：在第一课/第一课大作业/week1classification欢迎来到课程1的第一个作业！在这个任务中！您将通过使用Keras构建最先进的胸部X射线分类器来探索医学图像诊断。你将学会一下内容... 查看详情

第二课第一周1节-ai用于医学预后简介

第二门课程集中于医学预后（medicalprognosis）。预后是医学的一个分支，专门预测病人未来的健康状况。例如，根据病人的实验室结果，你能估计出未来5年内心脏病发作的风险吗？或是未来10年内死亡的风... 查看详情

第二课第一周3-4节-什么是预后？(代码片段)

什么是预后？我们首先要讨论什么是预后，以及为什么预后在医疗实践中很重要。预后是一个医学术语，指预测未来事件的风险。在这里，事件是一个通用术语，它描述了可能发生在个人身上的各种事情。事... 查看详情

第二课第一周2节-做医学预后，你需要掌握什么？

我想谈谈这门课的要求。本课程的设计重点在于为你提供概念和实用工具。你需要成功地为医学建立机器学习模型。对于这门课程，你不需要任何深入学习方法的背景知识，也不需要任何医学背景。不过，在你上这门... 查看详情

第二课第二周第1-5节-基于树的模型用于医学预后

本周，我们将使用决策树（DecisionTrees）构建我们的第一个机器学习模型。树(trees)在医学应用中非常有用的的原因是：1️⃣它们处理连续和分类数据的能力，2️⃣它们的可解释性以及训练速度。我们将使用树... 查看详情

吴恩达-医学图像人工智能专项课程-第一课第一周16-18节-如何确保数据集病人不重叠+作业解说(代码片段)

模型测试既然你已经了解了如何训练医学诊断模型，那么让我们来谈谈如何测试这样的模型。接下来你会学习如何测试这样的一个模型。您将学习如何正确使用训练、验证和测试集。以及为了评估你的模型需要强大的groundtrut... 查看详情

吴恩达-医学图像人工智能专项课程-第一课第一周6-10节总结+作业解读(代码片段)

...在医学图像分类问题上的一些前沿应用。本文将介绍第一课第一周6-10节的内容。主要讲解构建一个分类模型去识别胸片的肿块。以及分类模型将面临的三个挑战：类不平衡挑战、多任务挑战和数据集大小挑战。本节重点解决... 查看详情

收藏第一课第二周作业-学会计算分类各种指标-超详细教程(代码片段)

本次作业文件：在第一课/第一课大作业/week2metric这节课不需要对模型进行预测，所有的预测结果已经在csv文件中给出。作为提醒，我们的数据集包含14种不同情况的X射线，可通过X射线诊断。我们将使用我们在这... 查看详情

吴恩达-第一课第二周1-7节总结-医学深度学习模型的评估汇总(代码片段)

医学深度学习模型的评估汇总本周我们将深入探讨医学深度学习模型的评估。在医学上，由于决策具有很高的影响力，我们关心的是准确地了解模型何时对患者起作用，什么时候不起作用。您将学习一下指标，包... 查看详情