《GAIN》缺失值填补之源码解读

2,322 阅读8分钟

论文ICML 2018 《GAIN : Missing Data Imputation using Generative Adversarial Nets》 论文源码地址:

总览

这篇论文中作者提出了一种使用GAN来完成填补缺失值工作的GAIN,基本原理和标准的GAN相似,不同的在于根据具体问题所做的网络结构方面的改变。
GAIN中主要包括以下三个部分:

Generator,G:

它用来观察真实数据的每一部分,然后根据观察的结果填补缺失数据的部分,输出一个填补后完整的向量

Discriminator,D:

它接受一个完整的向量,来判别哪一部分数据是真实观察到的,哪一部分是被填补的

Hint Vector,hint:

它揭示了原始数据中缺失部分的某些信息,让D更加关注它所提示的部分,同时也逼迫G生成更加接近真实的数据用来填补。

一、参数配置:main_letter_spam.py 的 parse_args()

这个是我针对pycharm执行进行的修改,但是参数结构是一样的

def parse_args():
  # Inputs for the main function
  parser = argparse.ArgumentParser()
  parser.add_argument('--data_name',choices=['letter','spam'],default='spam',type=str)
  parser.add_argument('--miss_rate',help='missing data probability',default=0.2,type=float)
  parser.add_argument('--batch_size',help='the number of samples in mini-batch',default=128,type=int)
  parser.add_argument('--hint_rate',help='hint probability',default=0.9,type=float)
  parser.add_argument('--alpha',help='hyperparameter',default=100,type=float)
  parser.add_argument('--iterations',help='number of training interations',default=10000,type=int)
  return  parser.parse_args(args=[])

二、数据读取:datalloader.py

2.1数据集读取

if data_name in ['letter', 'spam']:
  file_name = 'data/'+data_name+'.csv'
  data_x = np.loadtxt(file_name, delimiter=",", skiprows=1)
elif data_name == 'mnist':
  (data_x, _), _ = mnist.load_data()
  data_x = np.reshape(np.asarray(data_x), [60000, 28*28]).astype(float)

这里的数据集只有letter和spam两个,如果不想下载数据集的话可以跑跑开源数据集mnist,

2.2数据集掩码(缺失)处理

# Parameters
no, dim = data_x.shape

# Introduce missing data
data_m = binary_sampler(1-miss_rate, no, dim)

miss_data_x = data_x.copy()

miss_data_x[data_m == 0] = np.nan

这里对数据集进行随机的掩码处理,可以理解为让他随机的把参数配置中的miss_rate指定值这里是20%的数数变为nan,也就是缺失值

  • 注:这里的binary_sampler是在untils中封装好的,很简单的一个随机选取
def binary_sampler(p, rows, cols):
  '''Sample binary random variables.
  
  Args:
    - p: probability of 1
    - rows: the number of rows
    - cols: the number of columns
    
  Returns:
    - binary_random_matrix: generated binary random matrix.
  '''
  unif_random_matrix = np.random.uniform(0., 1., size = [rows, cols])
  binary_random_matrix = 1*(unif_random_matrix < p)
  return binary_random_matrix

2.3返回数据集

return data_x, miss_data_x, data_m

返回数据的话,一共有三个:

返回数据数据含义
ori_data_x原始数据
miss_data_x缺失处理后数据
data_m掩码矩阵 其中0为缺失

三、gain模型

3.1数据归一化处理

# Normalization
norm_data, norm_parameters = normalization(data_x) #归一化

norm_data_x = np.nan_to_num(norm_data, 0) #nan值处理为0
  • data_x为原始数据
  • norm_data为归一化后的数据
  • norm_parameters为每一列的最大值数组,最小值数组这里叫做参数也没错,方便后续逆归一化

注:normalization也是untils中封装好的

def normalization (data, parameters=None):
  '''Normalize data in [0, 1] range.
  
  Args:
    - data: original data
  
  Returns:
    - norm_data: normalized data
    - norm_parameters: min_val, max_val for each feature for renormalization
  '''

  # Parameters
  _, dim = data.shape
  norm_data = data.copy()
  
  if parameters is None:
  
    # MixMax normalization
    min_val = np.zeros(dim)
    max_val = np.zeros(dim)
    
    # For each dimension
    for i in range(dim):
      min_val[i] = np.nanmin(norm_data[:,i])
      norm_data[:,i] = norm_data[:,i] - np.nanmin(norm_data[:,i])
      max_val[i] = np.nanmax(norm_data[:,i])
      norm_data[:,i] = norm_data[:,i] / (np.nanmax(norm_data[:,i]) + 1e-6)   
      
    # Return norm_parameters for renormalization
    norm_parameters = {'min_val': min_val,
                       'max_val': max_val}

  else:
    min_val = parameters['min_val']
    max_val = parameters['max_val']
    
    # For each dimension
    for i in range(dim):
      norm_data[:,i] = norm_data[:,i] - min_val[i]
      norm_data[:,i] = norm_data[:,i] / (max_val[i] + 1e-6)  
      
    norm_parameters = parameters    
      
  return norm_data, norm_parameters

3.2GAN网络基本结构

#  GAIN architecture
# Input placeholders
# Data vector
X = tf.placeholder(tf.float32, shape = [None, dim])
# Mask vector 
M = tf.placeholder(tf.float32, shape = [None, dim])
# Hint vector
H = tf.placeholder(tf.float32, shape = [None, dim])

一点前置准备,tf.placeholder()简单理解为占位符
为什么使用tf.placeholder(),就是因为可以在性能方面优化

Tensorflow的设计理念称之为计算流图,在编写程序时,首先构筑整个系统的graph,代码并不会直接生效,这一点和python的其他数值计算库(如Numpy等)不同,graph为静态的,类似于docker中的镜像。然后,在实际的运行时,启动一个session,程序才会真正的运行。这样做的好处就是:避免反复地切换底层程序实际运行的上下文,tensorflow帮你优化整个系统的代码。我们知道,很多python程序的底层为C语言或者其他语言,执行一行脚本,就要切换一次,是有成本的,tensorflow通过计算流图的方式,帮你优化整个session需要执行的代码,还是很有优势的。

所以placeholder()函数是在神经网络构建graph的时候在模型中的占位,此时并没有把要输入的数据传入模型,它只会分配必要的内存。等建立session,在会话中,运行模型的时候通过feed_dict()函数向占位符喂入数据。

3.2.1网络参数初始化

# 数据集的维度
no, dim = data_x.shape

# 隐藏单元维度
h_dim = int(dim)

# 判别器网络变量
D_W1 = tf.Variable(xavier_init([dim*2, h_dim])) # Data + Hint as inputs
D_b1 = tf.Variable(tf.zeros(shape = [h_dim]))
# print(D_W1.shape)
# print(D_W1)
D_W2 = tf.Variable(xavier_init([h_dim, h_dim]))
D_b2 = tf.Variable(tf.zeros(shape = [h_dim]))

D_W3 = tf.Variable(xavier_init([h_dim, dim]))
D_b3 = tf.Variable(tf.zeros(shape = [dim]))  # Multi-variate outputs

theta_D = [D_W1, D_W2, D_W3, D_b1, D_b2, D_b3]

#生成器网络变量
# Data + Mask as inputs (Random noise is in missing components)
G_W1 = tf.Variable(xavier_init([dim*2, h_dim]))  
G_b1 = tf.Variable(tf.zeros(shape = [h_dim]))

G_W2 = tf.Variable(xavier_init([h_dim, h_dim]))
G_b2 = tf.Variable(tf.zeros(shape = [h_dim]))

G_W3 = tf.Variable(xavier_init([h_dim, dim]))
G_b3 = tf.Variable(tf.zeros(shape = [dim]))

theta_G = [G_W1, G_W2, G_W3, G_b1, G_b2, G_b3]

3.2.1生成器generator

def generator(x,m):
  # 数据集和掩码矩阵拼接
  inputs = tf.concat(values = [x, m], axis = 1)

  G_h1 = tf.nn.relu(tf.matmul(inputs, G_W1) + G_b1)
  G_h2 = tf.nn.relu(tf.matmul(G_h1, G_W2) + G_b2)   
  # MinMax normalized output
  G_logit = tf.matmul(G_h2, G_W3) + G_b3
  G_prob = tf.nn.sigmoid(G_logit) 
  return G_prob
  • 注:tf.concat向量拼接axis的意义
# tensor t3 with shape [2, 3]
# tensor t4 with shape [2, 3]
tf.shape(tf.concat([t3, t4], axis = 0))  # [4, 3]
tf.shape(tf.concat([t3, t4], axis = 1))  # [2, 6]

3.2.2判别器discriminator

def discriminator(x, h):
  # 数据集和提示hint矩阵拼接
  inputs = tf.concat(values = [x, h], axis = 1) 
  D_h1 = tf.nn.relu(tf.matmul(inputs, D_W1) + D_b1)  
  D_h2 = tf.nn.relu(tf.matmul(D_h1, D_W2) + D_b2)
  D_logit = tf.matmul(D_h2, D_W3) + D_b3
  D_prob = tf.nn.sigmoid(D_logit)
  return D_prob
  • 这个generator生成器和discriminato判别器的结构完全一样,只有参数不同,至于为什么不同
  • hint揭示了原始数据中缺失部分的某些信息,其实和mask作用完全一样,但是discriminator判别器的hint的最终目的是逼近generator生成器的mask掩码矩阵,肯定不能直接把答案放到discriminator判别器网络里

3.3GAIN模型搭建

3.3.1生成器网络

## GAIN structure
# Generator 生成数据
G_sample = generator(X, M)

# Combine with observed data
Hat_X = X * M + G_sample * (1-M)

注:Hat_X为原始数据X的未缺失部分与生成网络生成的数据G_sample的缺失部分的合并

3.3.2判别器网络

# Discriminator
D_prob = discriminator(Hat_X, H)

3.3.3损失函数定义

## GAIN loss
D_loss_temp = -tf.reduce_mean(M * tf.log(D_prob + 1e-8) + (1-M) * tf.log(1. - D_prob + 1e-8)) 

G_loss_temp = -tf.reduce_mean((1-M) * tf.log(D_prob + 1e-8))

MSE_loss = tf.reduce_mean((M * X - M * G_sample)**2) / tf.reduce_mean(M)

D_loss = D_loss_temp
G_loss = G_loss_temp + alpha * MSE_loss 

D的损失函数就是判别器输出数据和M的差值的平均

  • 为什么不搞个绝对值,这样的差值有问题的吧?)
  • 经验证,没问题,sigmoid保证数据在0-1

G的损失函数包括两部分:(说实话一开始有点不太能理解,后续才悟到)

  • D生成的缺失部分数据的判定损失
  • 生成的未缺失值数据和真实的未缺失值数据的MSE的融合

注:

  • tf.log(y) 计算y的自然对数,如tf.log(1)=0
  • tf.reduce_mea函数用于计算张量tensor沿着指定的数轴(tensor的某一维度)上的平均值,主要用作降维或者计算tensor的平均值 如:tf.reduce_mean( input_tensor, axis=None, keep_dims=False, name=None, reduction_indices=None )
  • input_tensor: 输入的待降维的tensor
  • axis: 指定的轴,如果不指定,则计算所有元素的均值
  • keep_dims:是否保持维度,默认False。设置为True,输出的结果保持输入tensor的形状,设置为False,输出结果会降低维度
  • name: 操作的名称
  • reduction_indices:在以前版本中用来指定轴,已弃用

3.4模型训练

3.4.1初始化向量

## GAIN solver
D_solver = tf.train.AdamOptimizer().minimize(D_loss, var_list=theta_D)
G_solver = tf.train.AdamOptimizer().minimize(G_loss, var_list=theta_G)

## Iterations
sess = tf.Session()
sess.run(tf.global_variables_initializer())
  • 参数优化器AdamOptimizer,用来做学习率优化的
  • sess.run(tf.global_variables_initializer())在会话执行这个之后,之前的palceholder还有variable现在才初始化变量

3.4.2开始迭代训练

# Start Iterations
for it in tqdm(range(iterations)):    
    
  # Sample batch 小批次提取数据
  batch_idx = sample_batch_index(no, batch_size)
  X_mb = norm_data_x[batch_idx, :]  
  M_mb = data_m[batch_idx, :]  
  
  # Sample random vectors 
  Z_mb = uniform_sampler(0, 0.01, batch_size, dim) 
  
  # Sample hint vectors 
  H_mb_temp = binary_sampler(hint_rate, batch_size, dim)
  H_mb = M_mb * H_mb_temp
    
  # Combine random vectors with observed vectors
  X_mb = M_mb * X_mb + (1-M_mb) * Z_mb 
    
  _, D_loss_curr = sess.run([D_solver, D_loss_temp],feed_dict = {M: M_mb, X: X_mb, H: H_mb})
  _, G_loss_curr, MSE_loss_curr = sess.run([G_solver, G_loss_temp, MSE_loss],feed_dict = {X: X_mb, M: M_mb, H: H_mb})

注:

  • 这里的hint_rate照常理来说应该是1,但是这里用90%估计是个正则化设计,和dropout类似防过拟合的
  • 生成噪声数据Z_mb填补在缺失处
  • 利用训练好的生成器生成现在已经可以以假乱真的数据
  • 将生成的缺失部分的数据填补进去 sample_batch_index在unitls封装好了:
def sample_batch_index(total, batch_size):
  '''Sample index of the mini-batch.
  
  Args:
    - total: total number of samples
    - batch_size: batch size
    
  Returns:
    - batch_idx: batch index
  '''
  total_idx = np.random.permutation(total)
  batch_idx = total_idx[:batch_size]
  return batch_idx

uniform_sampler

def uniform_sampler(low, high, rows, cols):
  '''Sample uniform random variables.
  
  Args:
    - low: low limit
    - high: high limit
    - rows: the number of rows
    - cols: the number of columns
    
  Returns:
    - uniform_random_matrix: generated uniform random matrix.
  '''
  return np.random.uniform(low, high, size = [rows, cols])     

3.3.3 逆归一化

# Renormalization
imputed_data = renormalization(imputed_data, norm_parameters)  

# Rounding
imputed_data = rounding(imputed_data, data_x)  

renormalization

def renormalization (norm_data, norm_parameters):
  '''Renormalize data from [0, 1] range to the original range.
  
  Args:
    - norm_data: normalized data
    - norm_parameters: min_val, max_val for each feature for renormalization
  
  Returns:
    - renorm_data: renormalized original data
  '''
  
  min_val = norm_parameters['min_val']
  max_val = norm_parameters['max_val']

  _, dim = norm_data.shape
  renorm_data = norm_data.copy()
    
  for i in range(dim):
    renorm_data[:,i] = renorm_data[:,i] * (max_val[i] + 1e-6)   
    renorm_data[:,i] = renorm_data[:,i] + min_val[i]
    
  return renorm_data

rounding

def rounding (imputed_data, data_x):
  '''Round imputed data for categorical variables.
  
  Args:
    - imputed_data: imputed data
    - data_x: original data with missing values
    
  Returns:
    - rounded_data: rounded imputed data
  '''
  
  _, dim = data_x.shape
  rounded_data = imputed_data.copy()
  
  for i in range(dim):
    temp = data_x[~np.isnan(data_x[:, i]), i]
    # Only for the categorical variable
    if len(np.unique(temp)) < 20:
      rounded_data[:, i] = np.round(rounded_data[:, i])
      
  return rounded_data

最后经过数据处理返回即可,这个逆归一化很好理解,但是这个rounding我是真不好说这个unique后小于20的round成整数的实际含义。
输出了下,if len(np.unique(temp)) < 20:lette数据集没有执行的时候
这个可能是为了防止,全是缺失值的情况吧,毕竟我看到有很多log前+1e-06和+1e-08的

个人理解:

image.png

image.png

1.对于损失函数

D的损失函数就是判别器输出数据和M的差值的平均

D_loss_temp = -tf.reduce_mean(M * tf.log(D_prob + 1e-8) + (1-M) * tf.log(1. - D_prob + 1e-8))

G的损失函数包括两部分:

  • D生成的缺失部分数据的判定效果变差
  • 生成的未缺失值数据和真实的未缺失值数据的MSE变小

G_loss_temp = -tf.reduce_mean((1-M) * tf.log(D_prob + 1e-8))

MSE_loss = tf.reduce_mean((M * X - M * G_sample)**2) /tf.reduce_mean(M)

在迭代的后半阶段,D生成的缺失部分数据的判定损失一直在波动大多数时候维持在0.25左右,但是MSE一直在降

  1. 生成器G让生成的未缺失值数据逼近于真实数据,
  2. 然后如果能够保证判别器D不能识别出哪个是缺失的数据哪个是未缺失的数据,也就是说们的缺失值和未缺失值分布基本一致
  3. 至于我思考的为什么没有和缺失值数据的损失函数的问题,其实这个模型才是对实际场景的解决,也就是说我们平时的缺失数据都是真的缺失了不是模拟出来的,不能在训练的时候把真实数据搞进去(作者还是牛的,我陷入了思维误区)

2.对于网络架构

1. 生成器

  • 输入:数据矩阵X缺失标注矩阵M

缺失标注矩阵M:缺失处为0,反之为1

data_m = 1-np.isnan(data_x)

数据矩阵中X:未缺失部分是真实数据,缺失部分是噪声数据

X_mb = M_mb * X_mb + (1-M_mb) * Z_mb  

  • 输出:生成的数据样本,使未缺失的部分生成的数据逼近真实值

2. 鉴别器

  • 输入:数据矩阵 Hat_X缺失提示矩阵H

数据矩阵Hat_X:未缺失部分是真实数据,缺失部分是生成的数据

Hat_X = X * M + G_sample * (1-M)

缺失提示矩阵H:有一个有部分提示信息的标注矩阵,hint_rate=0.9

注:这是一种提示机制,强化了G和D的对抗过程。hint_rate不同那么给D的提示信息就不同,从而可以训练出具有多个分布的G,根据D的结果来选择最优的那一个。

H_mb_temp = binary_sampler(hint_rate, batch_size, dim)

随机把矩阵中90%的数据搞成1,剩下的为0,正则化设计?

H_mb = M_mb * H_mb_temp

  • 输出:对缺失值的判断,目标是逼近数据标注M