当前位置: 首页 > 编程日记 > 正文

Proximal Algorithms 4 Algorithms

目录

  • Proximal minimization
    • 解释
  • \(f(x) + g(x)\)
    • 解释1 最大最小算法
    • 不动点解释
    • Forward-backward 迭代解释
  • 加速 proximal gradient method
  • 交替方向方法 ADMM
    • 解释1 自动控制
    • 解释2 Augmented Largranians
    • 解释3 Flow interpretation
    • 解释4 不动点
    • 特别的情况 \(f(x) + g(Ax)\)

Proximal Algorithms

这一节介绍了一些利用proximal的算法.

Proximal minimization

这个相当的简单, 之前也提过,就是一个依赖不动点的迭代方法:
在这里插入图片描述
有些时候\(\lambda\)不是固定的:
\[ x^{k+1} := \mathbf{prox}_{\lambda^k f}(x^k), \sum_{k=1}^{\infty}\lambda^k = \infty \]

import numpy as np
import matplotlib.pyplot as plt 

\(f(x,y) = x^2 + 50y\)为例

f = lambda x: x[0] ** 2 + 50 * x[1] ** 2
x = np.linspace(-40, 40, 1000)
y = np.linspace(-20, 20, 500)
X, Y = np.meshgrid(x, y) #获取坐标
fig, ax = plt.subplots()
ax.contour(X, Y, f([X, Y]), colors="black")
plt.show()

在这里插入图片描述

求解proximal可得:
\[ x = \frac{v_1}{2\lambda + 1} \\ y = \frac{v_2}{100\lambda + 1} \]

def prox(v1, v2, lam):x = v1 / (2 * lam + 1)y = v2 / (100 * lam + 1)return x, y
times = 50
x = 30
y = 15
lam = 0.1
process = [(x, y)]
for i in range(times):x, y = prox(x, y, 0.1)process.append((x, y))
process = np.array(process)
x = np.linspace(-40, 40, 1000)
y = np.linspace(-20, 20, 500)
X, Y = np.meshgrid(x, y) #获取坐标
fig, ax = plt.subplots()
ax.contour(X, Y, f([X, Y]), colors="black")
ax.scatter(process[:, 0], process[:, 1])
ax.plot(process[:, 0], process[:, 1])
plt.show()

在这里插入图片描述

解释

除了之前已经提到过的一些解释:

Gradient flow

考虑下面的微分方程:
在这里插入图片描述
\(t \rightarrow \infty\)\(f(x(t))\rightarrow p^*\),其中\(p^*\)是最小值.
我们来看其离散的情形:
在这里插入图片描述

于是就有:
\[ x^{k+1} := x^k - h \nabla f(x^k) \]

还有一种后退的形式:
\[ \frac{x^{k+1}-x^k}{h}=-\nabla f(x^{k+1}) \]
此时,为了找到\(x^{k+1}\), 我们需要求解一个方程:
\[ x^{k+1} + h \nabla f(x^{k+1}) = x^k \\ \Rightarrow x^{k+1} = (I+ h \nabla f)^{-1}x^k = \mathbf{prox}_{hf}(x^k) \]

还有一种特殊的解释,这里不提了.

\(f(x) + g(x)\)

考虑下面的问题:
\[ \mathrm{minimize} \quad f(x) + g(x) \]
其中\(f\)是可微的.
我们可以通过下列proximal gradient method来求解:
\[ x^{k+1} := \mathbf{prox}_{\lambda^k g}(x^k - \lambda^k \nabla f(x^k)) \]
可以证明(虽然我不会),当\(\nabla f\) Lipschitz连续,常数为\(L\),那么,如果\(\lambda^k = \lambda \in (0, 1/L]\),这个方法会以\(O(1/k)\)的速度收敛.

还有一些直线搜素算法:
在这里插入图片描述
一般取\(\beta=1/2\)\(\widehat{f}_{\lambda}\)\(f\)的一个上界,在后面的解释中在具体探讨.

解释1 最大最小算法

最大最小算法, 最小化函数\(\varphi: \mathbb{R}^n \rightarrow \mathbb{R}\):
\[ x^{k+1} := \mathrm{argmin}_x \widehat{\varphi}(x, x^k) \]
其中\(\widehat{\varphi}(\cdot, x^k)\)\(\varphi\)的凸上界:\(\widehat{\varphi}(x, x^k) \ge \varphi(x)\), \(\widehat{\varphi}(x, x)=\varphi(x)\).
我们可以这么构造一个上界:
在这里插入图片描述
上面的式子很像泰勒二阶展开,首先这个函数符合第二个条件,下面我们证明,当\(\lambda \in (0, 1/L]\),那么它也符合第一个条件.
\[ \widehat{f}_{\lambda}(x) - f(x) = f(y) - f(x) +\nabla f(y)^T(x-y)+...=(\nabla f(y)-\nabla f(z))(x-y)+... \]
其中\(z = x + \theta (y-x), \theta \in [0, 1]\), 又Lipschitz连续,所以:
\[ \|\nabla f(y)-\nabla f(z)\| \le L\|y-z\|\le L\|y-x\| \]
考虑\(f(x+t\Delta x)\)关于\(t\)的二阶泰勒展式:
\[ f(x+t\Delta x) = f(x)+\nabla f(x)^T\Delta x t + \frac{1}{2}\Delta x^T \nabla^2f(x) \Delta x t^2 + o(t^2) \]
\(t=1\):
\[ f(x+\Delta x) = f(x)+\nabla f(x)^T\Delta x + \frac{1}{2}\Delta x^T \nabla^2f(x) \Delta x + ... \]
\[ \frac{\|\nabla f(x)-\nabla f(x+t\Delta x)\|}{t}\le L\|\Delta x\| \]
由当\(t \rightarrow 0\)时,左边为\(\|\nabla^2 f(x) \Delta x\|\), 所以\(\nabla^2 f(x)\)的最大特征值必小于\(L\), 所以:
\[ f(x+\Delta x) \le f(x)+\nabla f(x)^T\Delta x + \frac{L}{2}\|\Delta x\|_2^2 + ... \]
完蛋,好像只能证明在局部成立,能证明在全局成立吗?

\[ x^{k+1} := \mathrm{argmin}_x \widehat{f}_{\lambda}(x, x^k) \]
再令:
\[ q_{\lambda}(x, y)=\widehat{f}_{\lambda}(x,y) + g(x) \]
那么:
\[ x^{k+1} := \mathrm{argmin}_x q_{\lambda}(x, x^k)=\mathbf{prox}_{\lambda g}(x^k-\lambda \nabla f(x^k)) \]
上面的等式,可以利用第二节中的性质推出.

不动点解释

最小化\(f(x)+g(x)\)的点\(x^*\)应当满足:
\[ 0 \in \nabla f(x^*)+\partial g(x^*) \]
更一般地:
在这里插入图片描述
这便说明了一种迭代方式.

Forward-backward 迭代解释

考虑下列微分方程系统:
在这里插入图片描述
离散化后得:
在这里插入图片描述
注意,等式右边\(x^k\)\(x^{k+1}\),这正是巧妙之处.
解此方程可得:
在这里插入图片描述
这就是之前的那个迭代方法.

加速 proximal gradient method

其迭代方式为:
\[ y^{k+1} := x^k + w^k(x^k-x^{k-1}) \\ x^{k+1} := \mathbf{prox}_{\lambda^k g}(y^{k+1}-\lambda^k \nabla f(y^{k+1})) \]
\(w^k \in [0,1)\)
这个方法有点类似Momentum的感觉.
一个选择是:
\[ w^k = \frac{k}{k+3} \]
也有类似的直线搜索算法:

在这里插入图片描述

交替方向方法 ADMM

alternating direction method of multipliers (ADMM), 怎么说呢,久闻大名,不过还没看过类似的文章.
同样是考虑这个问题:
\[ \mathrm{minimize} \quad f(x) + g(x) \]
但是呢,这时\(f,g\)都不一定是可微的, ADMM采取的策略是:
\[ x^{k+1} := \mathbf{prox}_{\lambda f} (z^k - u^k) \\ z^{k+1} := \mathbf{prox}_{\lambda g} (x^{k+1} + u^k)\\ u^{k+1} := u^k + x^{k+1} -z^{k+1} \]

特殊的情况是, \(f\)\(g\)是指示函数,不妨设\(f\)是闭凸集\(\mathcal{C}\)的指示函数,而\(g\)是闭凸集\(\mathcal{D}\)的指示函数, 即:
\[ I_{\mathcal{C}}(x)=0, if \: x\in \mathcal{C}, else \: + \infty \]
这个时候,更新公式变为:
\[ x^{k+1} := \Pi_{\mathcal{C}} (z^k - u^k)\\ z^{k+1} := \Pi_{\mathcal{C}} (x^{k+1} + u^k) \\ u^{k+1} := u^k + x^{k+1} -z^{k+1} \]

解释1 自动控制

可以这么理解,\(z\)为状态,而\(u\)为控制,前俩步时离散时间动态系统(不懂啊...), 第三步的目标是选择\(u\)使得\(x=z\),所以\(x^{k+1}-z^{k+1}\)可以认为是一个信号误差,所以第三步就会把这些误差累计起来.

解释2 Augmented Largranians

我们可以将问题转化为:
在这里插入图片描述
augmented Largranian:
在这里插入图片描述
其中\(y\)为对偶变量.
\(z, y\)已知的条件下,最小化\(L\), 即:
\[ x^{k+1} := \mathrm{argmin}_x L_{\rho}(x, z^k, y^k) \]
\(x, y\)已知的条件下,最小化\(L\), 即:
\[ z^{k+1} := \mathrm{argmin}_z L_{\rho}(x^{k+1}, z, y^k) \]
最后一步:
\[ y^{k+1} := y^k + \rho (x^{k+1} - z^{k+1}) \]
如果依照对偶问题的知识,关于\(y\)应该是取最大,但是呢,关于\(y\)是一个仿射函数,所以没有最值,所以就简单地取那个?
注意到:
在这里插入图片描述
在这里插入图片描述
\(u^k = (1/\rho)y^k\), \(\lambda = 1/\rho\)就是最开始的结果.

解释3 Flow interpretation

问题(4.9)的最优条件(KKT条件):
在这里插入图片描述
其中\(v\)是对偶变量.考虑微分方程:
在这里插入图片描述
(4.11)取得稳定点的条件即为(4.10)(\(v=\rho)\)(这部分没怎么弄明白).
离散化情形为:
在这里插入图片描述
\(h = \lambda, \rho = 1/\lambda\)即可得ADMM.

解释4 不动点

原问题的最优条件为:
\[ 0 \in \partial f(x^*) + \partial g(x^*) \]

ADMM的不动点满足:
\[ x = \mathbf{prox}_{\lambda f} (x-u), \quad z = \mathbf{prox}_{\lambda g}(x+u), \quad u = u + x - z \]
从最后一个等式,我们可以知道:
\[ x = z \], 于是
\[ x = \mathbf{prox}_{\lambda f}(x - u), \quad x = \mathbf{prox}_{\lambda g}(x + u) \]
等价于:
\[ x = (I + \partial f)^{-1}(x - u), x = (I + \lambda \partial g)^{-1}(x + u) \]
等价于:
\[ x - u \in x + \lambda \partial f(x), \quad x + u \in x + \lambda \partial g(x) \]
俩个式子相加,说明\(x\)即为最优解.
再来说明一下,为什么可以相加,根据次梯度的定义:
\[ \lambda f(z) \ge \lambda f(x) + (-u)^T(z-x), \quad \forall z\in \mathbf{dom}f \\ \lambda g(z) \ge \lambda g(x) + (+u)^T(z-x), \quad \forall z\in \mathbf{dom}g \\ \]
相加可得:
\[ \lambda f(z) + \lambda g(z) \ge 2x + \lambda f(x) +\lambda g(x) + 0 \]
需要注意的是,我证明的时候也困扰了,
\[ x - u \in x + \lambda \partial f(x) \]
并不是指(x-u)是函数\(x^2/2 + \lambda f(x)\)的次梯度, 而是\(x-u\)\(\lambda f(x)\)的次梯度集合加上\(x\)的集合内,也就是\(-u\)是其次梯度.

对不起!又想当然了,其实没问题, 如果
\[ g \in \partial f_1(x) + h(x) \]
\(\partial f_2(x)=h(x)\)则:
\[ g \in \partial (f_1+f_2)(x) \]
证:
已知:
\[ f_1(z) \ge f_1(x)+\partial f_1(x)^T(z-x) \\ f_2(z) \ge f_2(x)+h(x)^T(z-x) \\ \]
俩式相加可得:
\[ (f_1+f_2)(z)\ge (f_1+f_2)(x) +(\partial f_1(x) + h(x))^T(z-x)=(f_1+f_2)(x) +g^T(z-x) \]
所以\(g \in \partial (f_1+f_2)(x)\), 注意\(g=g(x)\)也是无妨的.

特别的情况 \(f(x) + g(Ax)\)

考虑下面的问题:
\[ \mathrm{minimize} \quad f(x) + g(Ax) \]
上面的求解,也可以让\(\widetilde{g}(x) = g(Ax)\),这样子就可以用普通的ADMM来求解了, 但是有更加简便的方法.

在这里插入图片描述

这个的来源为:

在这里插入图片描述
再利用和之前一样的推导,不过,我要存疑的一点是最后的替代,我觉得应该是:
\[ \rho (A^TAx^k - A^T z^k)^T x + (1 / 2\mu) \|x-x^k\|_2^2 \]
否则推不出来啊.

转载于:https://www.cnblogs.com/MTandHJ/p/10994811.html

相关文章:

C# TripleDES NoPadding 时对待加密内容进行补字节(8个字节为一个Block)

补一个空格(半角): private static byte[] FormatData(String str) {var yu str.Length % 8;if (yu 0) return Encoding.GetEncoding(Consts.Charset).GetBytes(str);var size 8 - yu;var arr new byte[str.Length size];var data Enco…

keras Regressor 回归

import numpy as np np.random.seed(1337) # for reproducibility from keras.models import Sequential from keras.layers import Dense import matplotlib.pyplot as plt # 可视化模块import tensorflow as tf import keras.backend.tensorflow_backend as KTF# create som…

13、JsonResponse响应介绍

转载于:https://blog.51cto.com/yht1990/2406566

keras Classifier 分类

import numpy as np np.random.seed(1337) # for reproducibility from keras.models import Sequential from keras.layers import Dense, Activation from keras.optimizers import RMSprop import matplotlib.pyplot as plt # 可视化模块import tensorflow as tf import ke…

如何管理好自己的性格?

往往因为我们太感性,而获得与男人不一样的灵动的感受。而当过分的感性不合时宜地在职业生涯中表现出来时,我们该怎么调整自己呢?  由于女人与生俱来的特点,我们善良、有耐心,所以我们更易得到别人的支持和帮助&#…

Axis2 webservice入门--Webservice的发布与调用

一、Webservice发布 参考 http://www.cnblogs.com/demingblog/p/3263576.html 二、webservice 调用 部分参考:http://www.cnblogs.com/demingblog/p/3264688.html 使用myeclipse中的axis2插件生成客户端代码 new -->others到如下界面: 点next 到如下界…

Java断点续传(基于socket与RandomAccessFile的实现)

这是一个简单的C/S架构,基本实现思路是将服务器注册至某个空闲端口用来监视并处理每个客户端的传输请求。 客户端先获得用户给予的需传输文件与目标路径,之后根据该文件实例化RandomAccessFile为只读,之后客户端向服务器发送需传输的文件名文…

EJB调用原理分析

EJB调用原理分析 作者:robbin (MSN:robbin_fan AT hotmail DOT com) 版权声明:本文严禁转载,如有转载请求,请和作者联系 一个远程对象至少要包括4个class文件:远程对象;远程对象的接口;实现远程…

Jfinal Generator 不需要生成带某个前缀的表名数组的方法

2019独角兽企业重金招聘Python工程师标准>>> package com.demo.common.model; import javax.sql.DataSource; import com.jfinal.kit.PathKit; import com.jfinal.kit.Prop; import com.jfinal.kit.PropKit; import com.jfinal.plugin.activerecord.generato…

tensorflow 2

import tensorflow as tf import numpy as npdef test1():#create datax_datanp.random.rand(100).astype(np.float32)y_datax_data*0.10.3#create tensorflow structureWeightstf.Variable(tf.random_uniform([1],-1.0,1.0)) #一维,范围[-1,1]biasestf.Variable(tf…

PCB多层线路板打样难点

PCB多层板无论从设计上还是制造上来说,都比单双层板要复杂,一不小心就会遇到一些问题,那在PCB多层线路板打样中我们要规避哪些难点呢?  1、层间对准的难点  由于多层电路板中层数众多,用户对PCB层的校准要求越来越…

GARFIELD@11-07-2004

Vanity Fair转载于:https://www.cnblogs.com/rexhost/archive/2004/11/07/61286.html

python文件读写1

# -*- coding: utf-8 -*-# read txt file def readTextFile(file):f open(file, r)# 尽可能多的读取文件的内容,一般会将整个文件内容都会读取context f.read() print(context)f.close()def readTextFileByLines(file):f open(file, "r")lines f.read…

jfinal框架下使用c3P0连接池连接sql server 2008

2019独角兽企业重金招聘Python工程师标准>>> 闲话少说 进入正题 首先是工程需要的jar包 然后是c3p0的配置文件。我是这样配置的 仅供参考 jdbcDriver com.microsoft.sqlserver.jdbc.SQLServerDriver jdbcUrl jdbc:sqlserver://localhost:7777;databaseNametest us…

mongodb插入文档时不传ObjectId

type BookExt struct {ID bson.ObjectId bson:"_id"Title string bson:"title"SubTitle string bson:"subTitle"Author string bson:"author" } 以上结构体,在通过此结构体对象作为参数传入Insert插入…

[问题]DotNet 项目如何实现在构建时 Build 号自动增加?

[问题]DotNet 项目如何实现在构建时 Build 号自动增加? 继续昨天的问题,今天在Google上找了一下,没有找到很好的方案。目前找到的解决方案有以下几种:1.使用一个地三方的 VS.Net 插件,实现在编译时 Build 号自动增加&a…

编写程序记录文件位置

当我们编写程序是会注意到,首先是配置一些函数的结构体。 所以我们就要找到下面的界面,然后打开FWLB中.c文件下面所对应的.h文件,这样就能查找到相应的结构体。下图为我所找到的中断的结构体、 然后就是查找相对应的中断向量。具体就是打开 还…

mnist数据集保存为图片

#coding: utf-8 from tensorflow.examples.tutorials.mnist import input_data import scipy.misc import os import numpy as np# 读取MNIST数据集。如果不存在会事先下载。 mnist input_data.read_data_sets("MNIST_data/", one_hotTrue)# 我们把原始图片保存在MN…

Python3数据分析与挖掘建模实战

<div>课程地址&#xff1a;http://icourse8.com/Python3_shujufenxi.html</div>复制代码第1章 课程介绍【赠送相关电子书随堂代码】 第2章 数据获取 第3章 单因子探索分析与数据可视化 第4章 多因子探索分析 第5章 预处理理论 第6章 挖掘建模 第7章 模型评估 第8章…

tensorflow生成对抗网络

import tensorflow as tf import numpy as np import os from tensorflow.examples.tutorials.mnist import input_data from matplotlib import pyplot as pltBATCH_SIZE 64 UNITS_SIZE 128 LEARNING_RATE 0.001 EPOCH 300 SMOOTH 0.1print("mnist手写体生成对抗网络…

博客园今天早上是不是出现什么问题了?

下面是我进我的blog后台管理和浏览博客园给出的提示。大约几分钟后恢复正常。转载于:https://www.cnblogs.com/freeyzh/archive/2004/12/01/71269.html

模态框获取id一直不变,都是同一个id值

2019独角兽企业重金招聘Python工程师标准>>> $(.refund-btn).click(function(){//此处必须是$(this),否则$(.refund-btn)重新获取&#xff0c;导致值一直不变var id $(this).attr(data-id);//var id $(.refund-btn).attr(data-id);错误&#xff0c;这样会导致一直…

标准功能模块组件 -- 内部联络单组件,内部邮件组件,提高多人异地协同办公效率...

为什么80%的码农都做不了架构师&#xff1f;>>> 未必什么功能都需要自己开发&#xff0c;我们不会自己开发一个数据库系统&#xff0c;也不会自己开发一个操作系统&#xff0c;同样我们每个功能模块都未必需要自己开发&#xff0c;自己开发最核心的模块&#xff0c…

Microsoft patterns practices Enterprise Library released

一直关注这个东西&#xff0c;本来订阅了RSS&#xff0c;没想到GotDotNet上面的发布信息给清空了。 上周末发布的&#xff0c;今天才看到&#xff0c;刚刚下载了一个&#xff0c;下载还要求注册&#xff0c;真麻烦&#xff0c;现把地址共享&#xff0c;方便大家。 http://down…

图论之拓扑排序 poj 2367 Genealogical tree

题目链接 http://poj.org/problem?id2367 题意就是给定一系列关系&#xff0c;按这些关系拓扑排序。 #include<cstdio> #include<cstring> #include<queue> #include<vector> #include<algorithm> using namespace std; const int maxn200; int…

算法基础知识科普:8大搜索算法之顺序搜索

基本概念和术语 搜索表&#xff08;Search Table&#xff09;&#xff1a;是由同一类型的数据元素&#xff08;或记录&#xff09;构成的集合。 关键字&#xff08;Key&#xff09;&#xff1a;是数据元素中某个数据项的值&#xff0c;用它可以标识一个数据元素。若此关键字…

foj2024

为什么80%的码农都做不了架构师&#xff1f;>>> http://acm.fzu.edu.cn/problem.php?pid2024 View Code #include < stdio.h > #include < string .h > #define M 1010 int c[M][M]; int f[M][M]; int min( int a, int b, int c){ int z …

4701年新年快乐!

中华民族传统历法夏历&#xff08;农历&#xff09;采用的是干支纪年法&#xff0c;是世界上最古老的历法之一。干支即“六十甲子”&#xff0c;以60年为一循环。它的纪元开始相传可追溯到黄帝轩辕氏时代&#xff0c;按公元计算&#xff0c;第一个“甲子年”应是在公元前2697年…

Win10系列:JavaScript访问文件和文件夹

在实际开发中经常会遇到访问文件的情况&#xff0c;因此学习与文件有关的操作对程序开发很有帮助&#xff0c;关于文件操作的一些基本技术&#xff0c;在前面章节中有专门基于C#语言的详细讲解&#xff0c;本节主要介绍如何使用HTML5和JavaScript开发具有文件操作功能的Windows…

算法基础知识科普:8大搜索算法之二分搜索

昨天介绍了对无序搜素表的顺序搜索方法&#xff0c;今天介绍对有序搜索表的二分搜索方法&#xff0c;“二分”在算法设计中是非常常用的一种思想&#xff0c;除了处理如下普通的搜索外&#xff0c;还用于搜索方程的解等工程领域。但二分法仍然有缺陷&#xff0c;待后面慢慢介绍…