perturbed node penalty added (662f81e2) · Commits · 钟慕尧 / tvgl

ParallelTVGL.py

+15 −6

Original line number	Diff line number	Diff line
		@@ -73,6 +73,15 @@ def mp_parallel_tvgl((thetas, z0s, z1s, z2s, u0s, u1s, u2s,
		z0s[i] = pf.soft_threshold_odd(thetas[i] + u0s[i], lambd, rho)

		""" Z1-Z2 Update """
		if pen_func == "perturbed_node":
		for i in range(1, n):
		z1s[i-1], z2s[i] = pf.perturbed_node(thetas[i-1],
		thetas[i],
		u1s[i-1],
		u2s[i],
		beta,
		rho)
		else:
		for i in range(1, n):
		a = thetas[i] - thetas[i-1] + u2s[i] - u1s[i-1]
		e = getattr(pf, pen_func)(a, beta, rho)

SerialTVGL.py

+16 −10

Original line number	Diff line number	Diff line
		@@ -32,13 +32,19 @@ class SerialTVGL(TVGL):
		for i in range(self.blocks)]

		def z1_z2_update(self):
		if self.penalty_function == "perturbed_node":
		for i in range(1, self.blocks):
		self.z1s[i-1], self.z2s[i] = pf.perturbed_node(self.thetas[i-1],
		self.thetas[i],
		self.u1s[i-1],
		self.u2s[i],
		self.beta,
		self.rho)
		else:
		aa = [self.thetas[i] - self.thetas[i-1] + self.u2s[i] - self.u1s[i-1]
		for i in range(1, self.blocks)]
		ee = [getattr(pf, self.penalty_function)(a, self.beta, self.rho) for a in aa]
		#ee = pf.group_lasso_penaltys(aa, self.beta, self.rho)
		for i in range(1, self.blocks):
		# a = self.thetas[i] - self.thetas[i-1] + self.u2s[i] - self.u1s[i-1]
		# e = pf.group_lasso_penalty(a, 2*self.beta/self.rho)
		summ = self.thetas[i-1] + self.thetas[i] + self.u1s[i-1] + self.u2s[i]
		self.z1s[i-1] = 0.5*(summ - ee[i-1])
		self.z2s[i] = 0.5*(summ + ee[i-1])

penalty_functions.py

+53 −23

Original line number	Diff line number	Diff line
		@@ -6,8 +6,6 @@ def soft_threshold_odd(a, lambd, rho):
		parameter = lambd/rho
		dimension = np.shape(a)[0]
		e = np.eye(dimension)
		# = np.zeros((dimension, dimension))
		#p.fill_diagonal(e, np.diag(a))
		for i in range(dimension - 1):
		for j in range(i + 1, dimension):
		if abs(a[i, j]) > parameter:
		@@ -18,21 +16,6 @@ def soft_threshold_odd(a, lambd, rho):
		return e


		def soft_threshold_odds(aa, lambd, rho):
		parameter = lambd/rho
		dimension = np.shape(aa[0])[0]
		ee = [np.eye(dimension, dimension) for i in range(len(aa))]
		for a, e in zip(aa, ee):
		for i in range(dimension - 1):
		for j in range(i + 1, dimension):
		if abs(a[i, j]) > parameter:
		result = np.sign(a[i, j])*(
		abs(a[i, j]) - parameter)
		e[i, j] = result
		e[j, i] = result
		return ee


		def element_wise(a, beta, rho):
		nju = 2*beta/rho
		dimension = np.shape(a)[0]
		@@ -56,13 +39,60 @@ def group_lasso(a, beta, rho):
		return e


		def group_lasso_penaltys(aa, beta, rho):
		nju = 2*beta/rho
		dimension = np.shape(aa[0])[0]
		ee = [np.zeros((dimension, dimension)) for i in range(len(aa))]
		for a, e in zip(aa, ee):
		def perturbed_node(theta_pre, theta, u1, u2, beta, rho):

		""" Initialize ADMM algorithm """
		dimension = np.shape(theta)[0]
		nju = beta/(2*rho)
		y1 = np.ones((dimension, dimension))
		y2 = np.ones((dimension, dimension))
		v = np.ones((dimension, dimension))
		w = np.ones((dimension, dimension))
		uu1 = np.zeros((dimension, dimension))
		uu2 = np.zeros((dimension, dimension))
		c = np.zeros((dimension, 3*dimension))
		c[:, 0:dimension] = np.eye(dimension)
		c[:, dimension:2*dimension] = -np.eye(dimension)
		c[:, 2dimension:3dimension] = np.eye(dimension)
		ct = c.transpose()
		cc = np.linalg.inv(np.dot(ct, c) + 2np.eye(3dimension))

		""" Run algorithm """
		iteration = 0
		y_pre = []
		stopping_criteria = False
		while iteration < 1000 and stopping_criteria is False:
		a = (y1 - y2 - w - uu1 + (w.transpose() - uu2).transpose())/2

		""" V Update """
		e = np.zeros((dimension, dimension))
		for j in range(dimension):
		l2_norm = np.linalg.norm(a[:, j])
		if l2_norm > nju:
		e[:, j] = (1 - nju/l2_norm)*a[:, j]
		return ee
		v = e

		""" W, Y1, Y2 Update """
		b = np.zeros((3*dimension, dimension))
		b[0:dimension, :] = (v + uu2).transpose()
		b[dimension:2*dimension, :] = theta_pre + u1
		b[2dimension:3dimension, :] = theta + u2
		d = v + uu1
		wyy = np.dot(cc, 2*b - np.dot(ct, d))
		w = wyy[0:dimension, :]
		y1 = wyy[dimension:2*dimension, :]
		y2 = wyy[2dimension:3dimension, :]

		""" UU1, UU2 Update """
		uu1 = uu1 + (v + w) - (y1 - y2)
		uu2 = uu2 + v - w.transpose()

		""" Check stopping criteria """
		if iteration > 0:
		dif = y1 - y_pre
		fro_norm = np.linalg.norm(dif)
		if fro_norm < 1e-3:
		stopping_criteria = True
		y_pre = list(y1)
		iteration += 1
		return (y1, y2)

run.py

+3 −0

Original line number	Diff line number	Diff line
		@@ -14,6 +14,7 @@ if __name__ == "__main__" and len(sys.argv) > 1:
		# 2. Penalty function
		# 1 = "element_wise"
		# 2 = "group_lasso"
		# 3 = "perturbed_node
		# 3. Number of blocks to be created
		# 4. lambda
		# 5. beta
		@@ -34,6 +35,8 @@ if __name__ == "__main__" and len(sys.argv) > 1:
		penalty_function = "element_wise"
		elif sys.argv[2] == "2":
		penalty_function = "group_lasso"
		elif sys.argv[2] == "3":
		penalty_function = "perturbed_node"
		algo_type = "parallel"
		blocks = int(sys.argv[3])
		lambd = float(sys.argv[4])

Admin message