nlpfit.md

非线性拟合

先用

\[ f(x)=8\sin(x)\exp(x)-\frac{12}{\ln(x)} \]

生成数据。

然后给以上数据一些随机扰动：

\[ g(x) = f(x) + rnd \]

扰动\(rnd \in [-8000,8000]\)。

设置模型如下:

\[ p_{1}\sin(x)\exp(x)+\frac{p_{2}}{\ln(x)} \]

其解是\(p = [ 7.3725843942648,19.999999999999 ]^{T}\)，相关系数是\(0.97885363869532\)，残差平方和是\(464521063.49528\)。

随机扰动的取值范围越小，相关系数越接近于\(1\),\(p\)值也越接近\(8\)和\(-12\)，残差平方和越小。当扰动为零时，拟合模型与生成数据的\(f(x)\)完全重合。

local print = require('package.print'). print
local plt = require('package.plot')
local mol = libminoptlab

local xdat = mol.matrix.linspace(2,10,20)

function f(x)
	return 8*math.sin(x)*math.exp(x)-12/math.log(x)
end

local rnd = mol.random.new(-8000,8000)

function g(x)
	return f(x)+rnd:generate()
end

local ydat1 = mol.apply(xdat:tovector(),f)
local ydat2 = mol.apply(xdat:tovector(),g)

function model(x,p)
	return {math.sin(x[1])*math.exp(x[1]),1/math.log(x[1])}
end

local eqcon = mol.nullconfun()
local ineqcon = mol.nullconfun()

local lb = {0,-20}
local ub = {10,20}
local algo = mol.listoptalgo()[1]
local eps = 1e-5
local n = 1000

local p, corr,fval= mol.nlpfit(xdat,ydat2,model,lb,ub,eqcon,ineqcon,algo,eps,n)

local fitmodel = mol.tofit(model,p)
local fitx = mol.matrix.linspace(2,10,50)
local fity = fitx:apply(fitmodel)

print({
	p = p,
	corr = corr,
	fval = fval
})

local img1 = plt.scatter({{xdat:tovector(),ydat1}},{'f(x)'},{{marker='s',color='r'}})
local img2 = plt.scatter({{xdat:tovector(),ydat2}},{'g(x)'},{{marker='o',color='g'}})
local img3 = plt.line({{fitx:tovector(),fity}},{'fitmodel(x)'})
local fw = plt.framework(111,'非线性拟合',{'x','y'})
plt.show({{fw,{img1,img2,img3}}})