main_EYaleB_ori.py

import torch
import torch.nn as nn
import torch.nn.functional as F
import scipy.io as sio
from scipy.sparse.linalg import svds
from sklearn.preprocessing import normalize
from sklearn import metrics,cluster
import numpy as np
import random
from collections import OrderedDict
from utils import thrC,post_proC,Accuracy,purity,SelfExpression
import time
import math
import warnings

class Conv2dSamePad(nn.Module):
    """
    Implement Tensorflow's 'SAME' padding mode in Conv2d.
    When an odd number, say `m`, of pixels are need to pad, Tensorflow will pad one more column at right or one more
    row at bottom. But Pytorch will pad `m+1` pixels, i.e., Pytorch always pads in both sides.
    So we can pad the tensor in the way of Tensorflow before call the Conv2d module.
    """
    def __init__(self, kernel_size, stride):
        super(Conv2dSamePad, self).__init__()
        self.kernel_size = kernel_size if type(kernel_size) in [list, tuple] else [kernel_size, kernel_size]
        self.stride = stride if type(stride) in [list, tuple] else [stride, stride]

    def forward(self, x):
        in_height = x.size(2)
        in_width = x.size(3)
        out_height = math.ceil(float(in_height) / float(self.stride[0]))
        out_width = math.ceil(float(in_width) / float(self.stride[1]))
        pad_along_height = ((out_height - 1) * self.stride[0] + self.kernel_size[0] - in_height)
        pad_along_width = ((out_width - 1) * self.stride[1] + self.kernel_size[1] - in_width)
        pad_top = math.floor(pad_along_height / 2)
        pad_left = math.floor(pad_along_width / 2)
        pad_bottom = pad_along_height - pad_top
        pad_right = pad_along_width - pad_left
        return F.pad(x, [pad_left, pad_right, pad_top, pad_bottom], 'constant', 0)

class ConvTranspose2dSamePad(nn.Module):
    """
    This module implements the "SAME" padding mode for ConvTranspose2d as in Tensorflow.
    A tensor with width w_in, feed it to ConvTranspose2d(ci, co, kernel, stride), the width of output tensor T_nopad:
        w_nopad = (w_in - 1) * stride + kernel
    If we use padding, i.e., ConvTranspose2d(ci, co, kernel, stride, padding, output_padding), the width of T_pad:
        w_pad = (w_in - 1) * stride + kernel - (2*padding - output_padding) = w_nopad - (2*padding - output_padding)
    Yes, in ConvTranspose2d, more padding, the resulting tensor is smaller, i.e., the padding is actually deleting row/col.
    If `pad`=(2*padding - output_padding) is odd, Pytorch deletes more columns in the left, i.e., the first ceil(pad/2) and
    last `pad - ceil(pad/2)` columns of T_nopad are deleted to get T_pad.
    In contrast, Tensorflow deletes more columns in the right, i.e., the first floor(pad/2) and last `pad - floor(pad/2)`
    columns are deleted.
    For the height, Pytorch deletes more rows at top, while Tensorflow at bottom.
    In practice, we usually want `w_pad = w_in * stride` or `w_pad = w_in * stride - 1`, i.e., the "SAME" padding mode
    in Tensorflow. To determine the value of `w_pad`, we should pass it to this function.
    So the number of columns to delete:
        pad = 2*padding - output_padding = w_nopad - w_pad
    If pad is even, we can directly set padding=pad/2 and output_padding=0 in ConvTranspose2d.
    If pad is odd, we can use ConvTranspose2d to get T_nopad, and then delete `pad` rows/columns by
    ourselves.
    This module should be called after the ConvTranspose2d module with shared kernel_size and stride values.
    """

    def __init__(self, output_size):
        super(ConvTranspose2dSamePad, self).__init__()
        self.output_size = output_size

    def forward(self, x):
        in_height = x.size(2)
        in_width = x.size(3)
        pad_height = in_height - self.output_size[0]
        pad_width = in_width - self.output_size[1]
        pad_top = pad_height // 2
        pad_bottom = pad_height - pad_top
        pad_left = pad_width // 2
        pad_right = pad_width - pad_left
        return x[:, :, pad_top:in_height - pad_bottom, pad_left: in_width - pad_right]


class AE(nn.Module):
	def __init__(self,channel,kernel):
		super(AE, self).__init__()

		self.encoder = nn.Sequential(OrderedDict([
			('pad1',Conv2dSamePad(kernel[0],2)),
			('conv1',nn.Conv2d(channel[0],channel[1],kernel_size=kernel[0],stride=2)),
			('relu1',nn.ReLU()),

			('pad2',Conv2dSamePad(kernel[1],2)),
			('conv2',nn.Conv2d(channel[1],channel[2],kernel_size=kernel[1],stride=2)),
			('relu2',nn.ReLU()),

			('pad3',Conv2dSamePad(kernel[2],2)),
			('conv3',nn.Conv2d(channel[2],channel[3],kernel_size=kernel[2],stride=2)),
			('relu3',nn.ReLU())
			 ]))
		sizes = [[12, 11], [24, 21], [48, 42]]

		self.decoder = nn.Sequential(OrderedDict([
			('deconv1',nn.ConvTranspose2d(channel[3],channel[2],kernel_size=kernel[2],stride=2)),
			('padd1',ConvTranspose2dSamePad(sizes[0])),
			('relud1',nn.ReLU()),

			('deconv2',nn.ConvTranspose2d(channel[2],channel[1],kernel_size=kernel[1],stride=2)),
			('padd2',ConvTranspose2dSamePad(sizes[1])),
			('relud2',nn.ReLU()),

			('deconv3',nn.ConvTranspose2d(channel[1],channel[0],kernel_size=kernel[0],stride=2)),
			('padd3',ConvTranspose2dSamePad(sizes[2])),
			('relud3',nn.ReLU())
			 ]))

	def forward(self,x):
		z = self.encoder(x)
		output = self.decoder(z)     # shape=[n, c, w, h]

		return output

class EYaleB(nn.Module):
	def __init__(self,channel,kernel,num_sample):
		super(EYaleB,self).__init__()
		self.ae = AE(channel,kernel)
		self.n = num_sample
		self.self_expression = SelfExpression(self.n)
		self.self_expression_s = SelfExpression(self.n)

	def forward(self, x):
		z = self.ae.encoder(x)

		shape = z.shape
		z = z.view(self.n, -1)  # shape=[n, d]
		# pslb = self.fc(z)

		z_recon = self.self_expression(z)  # shape=[n, d]
		z_recon_reshape = z_recon.view(shape)

		x_recon = self.ae.decoder(z_recon_reshape)    # shape=[n, c, w, h]

		return x_recon, z, z_recon

	def loss_fn(self, x, x_recon, z, z_recon, weight_coef, weight_selfExp,weight_coef_s,weight_selfExp_s,w):
		cirtion=nn.MSELoss(reduction='sum')
		loss_ae = cirtion(x_recon, x)
		loss_coef = torch.norm(self.self_expression.Coefficient, p = 2) **2
		loss_coef_s = torch.norm(self.self_expression_s.Coefficient, p = 2) **2

		loss_selfExp = cirtion(z_recon, z)

		Coef = self.self_expression.Coefficient
		C = Coef
		loss_selfExp_s = cirtion(C, self.self_expression_s(C))

		contrastLoss = torch.zeros(1)
		Z1_matrix = F.normalize(z, dim=1)
		Y_matrix = F.normalize(z_recon,dim =1)
		temperature = 0.7

		cosin_similarity = (torch.matmul(Z1_matrix, Y_matrix.T)) 
		cosin_similarity = torch.exp(cosin_similarity / temperature)

		denominator = torch.sum(cosin_similarity,dim = 1)
		for i in range(denominator.shape[0]):
			contrastLoss += torch.log(cosin_similarity[i,i] / denominator[i])

		contrastLoss = w * contrastLoss       # 自监督

		# 总损失
		loss = loss_ae + weight_coef * loss_coef + weight_selfExp * loss_selfExp + weight_selfExp_s * loss_selfExp_s + weight_coef_s * loss_coef_s + contrastLoss

		return loss,loss_coef,loss_selfExp,loss_coef_s,loss_selfExp_s,contrastLoss

def train(model,x, y, epochs, lr=1e-3, weight_coef=1.0, weight_selfExp=150, weight_coef_s = 0.001,
	weight_selfExp_s = 1,w = 1,device='cuda',alpha=0.04, show=10):
    optimizer = torch.optim.Adam(model.parameters(), lr=lr)
    if not isinstance(x, torch.Tensor):
        x = torch.tensor(x, dtype=torch.float32, device=device)
    x = x.to(device)
    if isinstance(y, torch.Tensor):
        y = y.to('cpu').numpy()
    K = len(np.unique(y))

    log = "w:"+str(w)
    with open('./log/EYleB.txt','a') as f:
    	f.write('\n'+log)
    f.close()
    print(log)

    for epoch in range(epochs):
        train_tic = time.time()
        x_recon, z, z_recon = model(x)

        loss,loss_coef,loss_selfExp,loss_coef_s,loss_selfExp_s,contrastLoss = model.loss_fn(x, x_recon, z, z_recon,
         weight_coef=weight_coef,weight_selfExp=weight_selfExp,weight_coef_s = weight_coef_s,weight_selfExp_s=weight_selfExp_s,w = w)
        optimizer.zero_grad()
        loss.backward()
        optimizer.step()
        if epoch % show == 0 or epoch == epochs - 1:
            C1 = model.self_expression.Coefficient
            C2 = model.self_expression_s.Coefficient
            C = (C1 + C2).detach().to('cpu').numpy()
            Coef = thrC(C, alpha)
            tic = time.time()

            y_pred, C_post = post_proC(Coef, K, d, ro)

            toc = time.time()
            shijian = toc - tic
            print(shijian)

            acc = Accuracy(y, y_pred)
            pur = purity(y_pred,y)
            nmi = metrics.normalized_mutual_info_score(y, y_pred)
            ari = metrics.adjusted_rand_score(y, y_pred)


            train_info = 'Epoch %02d:loss=%.4f, acc = %.4f, nmi = %.4f, pur = %.4f, ari = %.4f' %(epoch,loss.item(), acc, nmi,pur,ari)
            train_loss = 'loss_coef = %.4f,loss_selfexp = %.4f,loss_coef_s = %.4f,loss_selfExp_s=%.4f,contrastLoss = %.4f'%(loss_coef.item(),loss_selfExp.item(),loss_coef_s.item(),loss_selfExp_s.item(),contrastLoss.item())
            with open('./log/EYleB.txt','a') as f:
            	f.write('\n'+train_info)
            	f.write('\n'+train_loss)
            f.close()
            print(train_info)
            print(train_loss)
        train_toc = time.time()
        print(train_toc - train_tic)

if __name__ == "__main__":
	# load data
	data = sio.loadmat('Data/YaleBCrop025.mat')
	img = data['Y']
	I = []
	Label = []
	for i in range(img.shape[2]):
	    for j in range(img.shape[1]):
	        temp = np.reshape(img[:, j, i], [42, 48])
	        Label.append(i)
	        I.append(temp)
	I = np.array(I)
	y_total = np.array(Label[:])
	Img = np.transpose(I, [0, 2, 1])
	x_total = np.expand_dims(Img[:], 1)
	print(y_total)

	num_class = 38

	# network and optimization parameters
	num_sample = num_class * 64
	channel = [1, 10, 20, 30]
	kernel = [5,3,3]
	num_epoch = 1500
	weight_coef = 1
	weight_selfExp = 1.0 * 10 ** (num_class / 10.0 - 3.0)
	weight_coef_s = 1
	weight_selfExp_s = 0.01
	w = [0.001,0.01,0.1,1,10,100,1000]
	d = 11
	ro = 2.5

	# post clustering parameters
	alpha = 0.09  # threshold of C
	device = torch.device('cuda' if torch.cuda.is_available() else 'cpu')

	random.seed(1)
	seed = np.random.seed(1)
	torch.manual_seed(1)
	torch.backends.cudnn.benchmark = False
	torch.backends.cudnn.deterministic = True

	for i in w:
		model=EYaleB(channel,kernel,num_sample)
		model.to(device)
		ae_weight = torch.load("./pretrained-EYaleB/yaleb.pkl")
		model.ae.load_state_dict(ae_weight)
		print("Pretrained ae weights are loaded successfully.")
		warnings.filterwarnings("ignore")
		train(model, x_total, y_total, num_epoch, weight_coef=weight_coef, weight_selfExp=weight_selfExp,
			weight_coef_s = weight_coef_s,weight_selfExp_s = weight_selfExp_s,w = i,alpha=alpha, show=50, device=device)