zoukankan      html  css  js  c++  java
  • 优化算法——拟牛顿法之L-BFGS算法

    一、BFGS算法

        在“优化算法——拟牛顿法之BFGS算法”中,我们得到了BFGS算法的校正公式:




    利用Sherman-Morrison公式可对上式进行变换,得到




    ,则得到:




    二、BGFS算法存在的问题

        在BFGS算法中。每次都要存储近似Hesse矩阵,在高维数据时,存储浪费非常多的存储空间,而在实际的运算过程中。我们须要的是搜索方向。因此出现了L-BFGS算法。是对BFGS算法的一种改进算法。

    L-BFGS算法中。仅仅保存近期的次迭代信息。以减少数据的存储空间。

    三、L-BFGS算法思路

        令,则BFGS算法中的能够表示为:




    若在初始时,假定初始的矩阵,则我们能够得到:









    若此时。仅仅保留近期的步:




    这样在L-BFGS算法中。不再保存完整的。而是存储向量序列。须要矩阵时,使用向量序列计算就能够得到。而向量序列也不是全部都要保存,仅仅要保存最新的步向量就可以。

    四、L-BFGS算法中的方向的计算方法


    五、实验仿真

    lbfgs.py

    #coding:UTF-8
    
    from numpy import *
    from function import *
    
    def lbfgs(fun, gfun, x0):
        result = []#保留终于的结果
        maxk = 500#最大的迭代次数
        rho = 0.55
        sigma = 0.4
        
        H0 = eye(shape(x0)[0])
        
        #s和y用于保存近期m个,这里m取6
        s = []
        y = []
        m = 6
        
        k = 1
        gk = mat(gfun(x0))#计算梯度
        dk = -H0 * gk
        while (k < maxk):             
            n = 0
            mk = 0
            gk = mat(gfun(x0))#计算梯度
            while (n < 20):
                newf = fun(x0 + rho ** n * dk)
                oldf = fun(x0)
                if (newf < oldf + sigma * (rho ** n) * (gk.T * dk)[0, 0]):
                    mk = n
                    break
                n = n + 1
            
            #LBFGS校正
            x = x0 + rho ** mk * dk
            #print x
            
            #保留m个
            if k > m:
                s.pop(0)
                y.pop(0)
                
            #计算最新的
            sk = x - x0
            yk = gfun(x) - gk
            
            s.append(sk)
            y.append(yk)
            
            #two-loop的过程
            t = len(s)
            qk = gfun(x)
            a = []
            for i in xrange(t):
                alpha = (s[t - i - 1].T * qk) / (y[t - i - 1].T * s[t - i - 1])
                qk = qk - alpha[0, 0] * y[t - i - 1]
                a.append(alpha[0, 0])
            r = H0 * qk
                
            for i in xrange(t):
                beta = (y[i].T * r) / (y[i].T * s[i])
                r = r + s[i] * (a[t - i - 1] - beta[0, 0])
    
                
            if (yk.T * sk > 0):
                dk = -r            
            
            k = k + 1
            x0 = x
            result.append(fun(x0))
        
        return result
    

    function.py

    #coding:UTF-8
    '''
    Created on 2015年5月19日
    
    @author: zhaozhiyong
    '''
    
    from numpy import *
    
    #fun
    def fun(x):
        return 100 * (x[0,0] ** 2 - x[1,0]) ** 2 + (x[0,0] - 1) ** 2
    
    #gfun
    def gfun(x):
        result = zeros((2, 1))
        result[0, 0] = 400 * x[0,0] * (x[0,0] ** 2 - x[1,0]) + 2 * (x[0,0] - 1)
        result[1, 0] = -200 * (x[0,0] ** 2 - x[1,0])
        return result
    

    testLBFGS.py

    #coding:UTF-8
    '''
    Created on 2015年6月6日
    
    @author: zhaozhiyong
    '''
    
    from lbfgs import *
    
    import matplotlib.pyplot as plt  
    
    x0 = mat([[-1.2], [1]])
    result = lbfgs(fun, gfun, x0)
    print result
    
    n = len(result)
    ax = plt.figure().add_subplot(111)
    x = arange(0, n, 1)
    y = result
    ax.plot(x,y)
    
    plt.show()

    实验结果



    參考文献


  • 相关阅读:
    linuxepoll研究 Geek_Ma 博客园
    socklen_t 类型 blueliuyun的专栏 博客频道 CSDN.NET
    自己动手写web服务器一(浏览器的访问信息) 任天胜的个人空间 开源中国社区
    UNIX Domain Socket IPC blueliuyun的专栏 博客频道 CSDN.NET
    How to use epoll? A complete example in C
    gzip头部格式 任天胜的个人空间 开源中国社区
    CWnd与HWND的区别与转换
    MFC 框架各部分指针获取方式
    windows 注册表的编程
    VS2010创建C++项目类向导和智能感知不可用
  • 原文地址:https://www.cnblogs.com/zhchoutai/p/6803832.html
Copyright © 2011-2022 走看看