飞桨黑客马拉松第四期 科学计算 Learning to regularize with a variational autoencoder for hydrologic inverse analysis

1.简介

本项目基于paddle框架复现。

论文主要点如下：

• 作者提出了一种基于变分自动编码器 (VAE)的正则化方法；
• 这种方法的优点1: 对来自VAE的潜在变量（此处指encoder的输出）执行正则化，可以简单地对其进行正则化；
• 这种方法的优点2: VAE减少了优化问题中的变量数量，从而在伴随方法不可用时使基于梯度的优化在计算上更加高效。

本项目关键技术要点：

• 实现paddle和julia混写代码梯度传递，避免大面积重写julia代码并能够调用优秀的julia代码；
• 发现paddle minimize_lbfgs存在问题, 待提交issue确认。

实验结果要点：

• 成功复现论文代码框架及全流程运行测试；
• 本次复现精度因无法使用相同样本，无法与论文中数据进行一一比较。本项目给出了采用paddle编写的框架结果。

论文信息： O'Malley D, Golden J K, Vesselinov V V. Learning to regularize with a variational autoencoder for hydrologic inverse analysis[J]. arXiv preprint arXiv:1906.02401, 2019.

参考GitHub地址： https://github.com/madsjulia/RegAE.jl

项目aistudio地址： https://aistudio.baidu.com/aistudio/projectdetail/5541961

模型结构

2.数据集

本项目数据集通过作者提供的julia代码生成，生成后保存为npz文件，已上传aistudio数据集并关联本项目。

以下为数据生成和数据保存代码的说明

（1）作者通过julia中的DPFEHM和GaussianRandomFields进行数据生成，代码可参考本项目/home/aistudio/RegAE.jl/examples/hydrology/ex_gaussian.jl，可根据其中参数进行修改；

（2）数据保存代码。在/home/aistudio/RegAE.jl/examples/hydrology/ex.jl代码中增加以下代码，可将数据通过转换为numpy数据并保存为npz。

using Distributed
using PyCall # 增加此处引用

@everywhere variablename = "allloghycos"
@everywhere datafilename = "$(results_dir)/trainingdata.jld2"
if !isfile(datafilename)
    if nworkers() == 1
        error("Please run in parallel: julia -p 32")
    end
    numsamples = 10^5
    @time allloghycos = SharedArrays.SharedArray{Float32}(numsamples, ns[2], ns[1]; init=A -> samplehyco!(A; setseed=true))
    # @time JLD2.@save datafilename allloghycos
    
    ########### 此处为增加部分 ###########
    p_trues = SharedArrays.SharedArray{Float32}(3, ns[2], ns[1]; init=samplehyco!) # 计算p_true

    np = pyimport("numpy")
    training_data = np.asarray(allloghycos)
    test_data = np.asarray(p_trues)

    np_coords = np.asarray(coords)
    np_neighbors = np.asarray(neighbors)
    np_areasoverlengths = np.asarray(areasoverlengths)
    np_dirichletnodes = np.asarray(dirichletnodes)
    np_dirichletheads = np.asarray(dirichletheads)

    np.savez("$(results_dir)/gaussian_train.npz",
        data=training_data,
        test_data=test_data,
        coords=np_coords,
        neighbors=np_neighbors,
        areasoverlengths=np_areasoverlengths,
        dirichletnodes=np_dirichletnodes,
        dirichletheads=np_dirichletheads)
end

3.环境依赖

本项目为julia和python混合项目。

julia依赖

• DPFEHM
• Zygote

python依赖

• paddle
• julia (pip安装)
• matplotlib

本项目已经提供安装后压缩文档，可fork本项目后执行以下代码进行解压安装。

# 解压预下载文件和预编译文件
!tar zxf /home/aistudio/opt/curl-7.88.1.tar.gz -C /home/aistudio/opt # curl 预下载文件
!tar zxf /home/aistudio/opt/curl-7.88.1-build.tgz -C /home/aistudio/opt # curl 预编译文件
!tar zxf /home/aistudio/opt/julia-1.8.5-linux-x86_64.tar.gz -C /home/aistudio/opt # julia 预下载文件
!tar zxf /home/aistudio/opt/julia_package.tgz -C /home/aistudio/opt # julia依赖库文件
!tar zxf /home/aistudio/opt/external-libraries.tgz -C /home/aistudio/opt # pip依赖库文件

####### 以下指令需要时可参考执行，上述压缩包已经完成以下内容 #######

# curl 编译指令，当解压后无效使用
!mkdir -p /home/aistudio/opt/curl-7.88.1-build
!/home/aistudio/opt/curl-7.88.1/configure --prefix=/home/aistudio/opt/curl-7.88.1-build --with-ssl --enable-tls-srp
!make install -j4 

# 指定curl预编译文件
!export LD_LIBRARY_PATH=/home/aistudio/opt/curl-7.88.1-build/lib:$LD_LIBRARY_PATH
!export PATH=/home/aistudio/opt/curl-7.88.1-build/bin:$PATH
!export CPATH=/home/aistudio/opt/curl-7.88.1-build/include:$CPATH
!export LIBRARY_PATH=/home/aistudio/opt/curl-7.88.1-build/lib:$LIBRARY_PATH

# 指定已经安装的julia包
!export JULIA_DEPOT_PATH=/home/aistudio/opt/julia_package
# 指定julia使用清华源
!export JULIA_PKG_SERVER=https://mirrors.tuna.tsinghua.edu.cn/julia
# julia 安装依赖库
# 需要先export JULIA_DEPOT_PATH 环境变量，否则安装位置为~/.julia, aistudio无法保存
!/home/aistudio/opt/julia-1.8.5/bin/julia -e "using Pkg; Pkg.add(\"DPFEHM\")" 
!/home/aistudio/opt/julia-1.8.5/bin/julia -e "using Pkg; Pkg.add(\"Zygote\")"
!/home/aistudio/opt/julia-1.8.5/bin/julia -e "using Pkg; Pkg.add(\"PyCall\")"

使用方法可以参考以下代码和julia导数传递测试.ipynb文件。

import paddle
import os 
import sys

# julia 依赖
os.environ['JULIA_DEPOT_PATH'] = '/home/aistudio/opt/julia_package'
# pip 依赖
sys.path.append('/home/aistudio/opt/external-libraries')

# julieries
from julia.api import Julia

jl = Julia(compiled_modules=False,runtime="/home/aistudio/opt/julia-1.8.5/bin/julia")
# import julia
from julia import Main

4.快速开始

本项目运行分为两个步骤：

• （1）训练步骤。通过运行train.ipynb文件，可以得到训练后的模型参数，具体代码请参考train.ipynb文件及其中注释说明；
• （2）测试步骤。通过运行test.ipynb文件，应用训练后的模型参数，对latent domain进行优化。

5.代码结构与详细说明

├── data                        #预生成数据文件
│   └── data193595
├── main.ipynb                  #本说明文件
├── opt                         #环境配置文件，已压缩，解压即可使用
│   ├── curl-7.88.1
│   ├── curl-7.88.1-build
│   ├── curl-7.88.1-build.tgz
│   ├── curl-7.88.1.tar.gz
│   ├── external-libraries
│   ├── external-libraries.tgz
│   ├── julia-1.8.5
│   ├── julia-1.8.5-linux-x86_64.tar.gz
│   ├── julia_package
│   └── julia_package.tgz
├── params_vae_nz100            #模型参数文件
│   └── model.pdparams
├── params_vae_nz200
│   └── model.pdparams
├── params_vae_nz400
│   └── model.pdparams
├── test.ipynb                  #测试文件
├── train.ipynb                 #训练文件
├── julia导数传递测试.ipynb       #julia和python混合测试文件

train文件和test文件关联性说明

我们依照论文作者的符号进行说明，$p$为数据输入，$\hat{p}$为数据输入，$loss=mse(p,\hat{p}) + loss_{kl}(\hat{p},N(0,1))$。

• （1）通过train能够得到训练后的Autoencoder(包含encoder和decoder)；
• （2）通过test调用训练后的encoder针对testdata生成latent_test，并得到latent_mean；
• （3）针对新生成的样本$p_{new}$，通过LBFGS方法不断优化latent_mean，直到obj_fun最小，其中obj_fun = mse($p_{new}$,$\hat{p}{new}$)+mse(sci_fun($p{new}$),sci_fun($\hat{p}_{new}$))，sci_fun为任何其他科学计算模拟方法。

paddle.incubate.optimizer.functional.minimize_lbfgs 问题

以下为paddle官方minimize_lbfgs API:

paddle.incubate.optimizer.functional.minimize_lbfgs(objective_func, initial_position, history_size=100, max_iters=50, tolerance_grad=1e-08, tolerance_change=1e-08, initial_inverse_hessian_estimate=None, line_search_fn='strong_wolfe', max_line_search_iters=50, initial_step_length=1.0, dtype='float32', name=None)

• （1）参数max_line_search_iters无效。虽然设置了此参数，但是内部没有传递对应参数；
• （2）中wolfe条件1错误。line256处应为phi_2 >= phi_1，以下为paddle部分源码。

        # 1. If phi(a2) > phi(0) + c_1 * a2 * phi'(0) or [phi(a2) >= phi(a1) and i > 1], 
        #         a_star= zoom(a1, a2) and stop;
        pred1 = ~done & ((phi_2 > phi_0 + c1 * a2 * derphi_0) |
                         ((phi_2 >= phi_0) & (i > 1)))

6.复现结果

不同latent维度对比

通过实验结果可以发现：

• （1）不同样本之间存在差距，并不是所有样本都能优化得到良好的latent变量；
• （2）随着模型latent维度的上升，模型效果逐渐提升。

latent_random和latent_mean对比

本项目还增加了latent_random和latent_mean对生成结果的对比。此处对latent_random和latent_mean再次说明：

• latent_random：通过paddle.randn生成的高斯噪声得到；
• latent_mean：通过对所有testdata进行encoder结果平均得到。

以下为通过latent_random得到的实验结果

通过对比，可以发现latent_mean对优化结果重要影响。近似正确的latent变量能够得到更优的生成结果。

LBFGS优化收敛情况

可以从如下图中看出，使用paddle minimize_lbfgs能够有效优化收敛。

7.延伸思考

如果深入思考本项目，会发现模型在test过程中是使用真实数据作为目标进行lbfgs优化，这种计算方式还有意义吗？

回答是肯定的！有意义！

以下为本人个人观点：

• （1）通过实验对比latent_random和latent_mean的最终生成结果差距，可以发现一个良好的初值对模型的影响是巨大的。当前diffusion模型在sample生成过程对生成的高斯噪声不断做denoise操作，这其中生成的噪声数据如果经过预先优化，不仅能够加速diffusion的生成速度，而且能够提升数据的生成质量。
• （2）在域迁移等研究领域，可以使用这种latent逐渐生成中间过渡变量，达到不同域数据的迁移生成。

7.模型信息

信息	说明
发布者	朱卫国 (DrownFish19)
发布时间	2023.03
框架版本	paddle 2.4.1
支持硬件	GPU、CPU
aistudio	notebook

请点击此处查看本环境基本用法.
Please click here for more detailed instructions.