 |
 |
 |
 |
 |
 |
||
 |
 |
 |
 |
 |
| 3 | 1/1 | 返回列表 |
| 查看: 845 | 回復(fù): 2 | ||
manmanbobo新蟲 (初入文壇)
|
[求助]
LM算法碰到分段函數(shù)該怎么處理 已有1人參與
|
|
各位蟲友大家好，本蟲在目前在用LM算法進(jìn)行函數(shù)的曲線擬合，函數(shù)分兩段，一部分類似于余弦，另一部分類似于冪函數(shù)。在循環(huán)中，因為分段點(diǎn)也是一個未知的變量，所以不能進(jìn)行比較大小，特地求助廣大蟲友，感激不盡。這個問題卡了我好久。謝謝。下面附上一些代碼。 代碼1：單一函數(shù)的擬合，為了蟲友能更好理解LM算法。 % 計算函數(shù)f的雅克比矩陣，是解析式 syms a b c d y x real; f=a+b*cos((pi/c)*(x-d)); Jsym=jacobian(f,[a b c d]) % 擬合用數(shù)據(jù)。 data_1=[0.25 0.5 1 1.5 2 3 4 6 8]; obs_1=[19.21 18.15 15.36 14.10 12.89 9.32 7.45 5.24 3.01]; % 2. LM算法 % 初始猜測s a0=4; b0=20;c0=14;d0=13; y_init = a0+b0*cos((pi/c0)*(x-d0)); % 數(shù)據(jù)個數(shù) Ndata=length(obs_1); % 參數(shù)維數(shù) Nparams=4; % 迭代最大次數(shù) n_iters=50; % LM算法的阻尼系數(shù)初值 lamda=0.01; % step1: 變量賦值 updateJ=1; a_est=a0; b_est=b0; c_est=c0; d_est=d0; % step2: 迭代 for it=1:n_iters if updateJ==1 % 根據(jù)當(dāng)前估計值，計算雅克比矩陣 J=zeros(Ndata,Nparams); J for i=1:length(data_1) i J(i,  =[ 1 cos((pi*(d_est - data_1))/c_est) (b_est*pi*sin((pi*(d_est - data_1))/c_est).*(d_est - data_1))/c_est^2 -(b_est*pi*sin((pi*(d_est - data_1))/c_est))/c_est];end % 根據(jù)當(dāng)前參數(shù)，得到函數(shù)值 y_est =a_est+b_est*cos((pi/c_est)*(data_1-d_est)) ; % 計算誤差 d=obs_1-y_est; % 計算（擬）海塞矩陣 H=J'*J; % 若是第一次迭代，計算誤差 if it==1 e=dot(d,d); end end % 根據(jù)阻尼系數(shù)lamda混合得到H矩陣 H_lm=H+(lamda*eye(Nparams,Nparams)); % 計算步長dp，并根據(jù)步長計算新的可能的\參數(shù)估計值 dp=inv(H_lm)*(J'*d( );g = J'*d( ;a_lm=a_est+dp(1); b_lm=b_est+dp(2); c_lm=c_est+dp(3); d_lm=d_est+dp(4); % 計算新的可能估計值對應(yīng)的y和計算殘差e y_est_lm =a_lm+b_lm*cos((pi/c_lm)*(data_1-d_lm)) ; d_lm=obs_1-y_est_lm; e_lm=dot(d_lm,d_lm); % 根據(jù)誤差，決定如何更新參數(shù)和阻尼系數(shù) if e_lm<e lamda=lamda/10; a_est=a_lm; b_est=b_lm; c_est=c_lm; d_est=d_lm; e=e_lm; disp(e); updateJ=1; else updateJ=0; lamda=lamda*10; end end %顯示優(yōu)化的結(jié)果 a_est b_est c_est d_est 代碼2（部分）：對代碼1進(jìn)行改進(jìn)，想求分段函數(shù)的偏導(dǎo)。 % 計算函數(shù)f的雅克比矩陣，是解析式 syms a b c d f g x k y1 y2 real; if x<f y1=a+b*cos((pi/c)*(x-d)); y1x=jacobian(y1,[a b c d]); else y2=(a+g)+[b*cos((pi/c)*(f-d))-g]*exp(-(x-f)/k); y2x=jacobian(y2,[a b c d g ]) end %Jsym=jacobian(f,[a b c d f g h]); % 擬合用數(shù)據(jù)。 data_1=[0.25 0.5 1 1.5 2 3 4 6 8]; obs_1=[19.21 18.15 15.36 14.10 12.89 9.32 7.45 5.24 3.01]; % 2. LM算法 % 初始猜測s a0=4; b0=20;c0=14;d0=13;f0=17;g0=0.9; y_init = a0+b0*cos((pi/c0)*(x-d0)); % 數(shù)據(jù)個數(shù) Ndata=length(obs_1); % 參數(shù)維數(shù) Nparams=4; % 迭代最大次數(shù) n_iters=50; % LM算法的阻尼系數(shù)初值 lamda=0.01; % step1: 變量賦值 updateJ=1; a_est=a0; b_est=b0; c_est=c0; d_est=d0; f_est=f0; g_est=g0; % step2: 迭代 for it=1:n_iters if updateJ==1 % 根據(jù)當(dāng)前估計值，計算雅克比矩陣 J=zeros(Ndata,Nparams); a=ones(9,1); if x<h J=[ a cos((pi*(d_est - data_1))/c_est)' ((b_est*pi*sin((pi*(d_est - data_1))/c_est).*(d_est - data_1))/c_est^2)' (-(b_est*pi*sin((pi*(d_est - data_1))/c_est))/c_est)'] else J=[ a cos((pi*(d - f))/c)*exp((f - x)/k), (b*pi*exp((f - x)/k)*sin((pi*(d - f))/c)*(d - f))/c^2, -(b*pi*exp((f - x)/k)*sin((pi*(d - f))/c))/c, 1 - exp((f - x)/k), (b*pi*exp((f - x)/k)*sin((pi*(d - f))/c))/c] end 再次感謝蟲友，謝謝！ |
回復(fù)此樓» 猜你喜歡
2026年循環(huán)經(jīng)濟(jì)功能材料國際會議（ICFMCE 2026）
已經(jīng)有0人回復(fù)
-
2026年第五屆電氣、電子與信息工程國際會議（ISEEIE 2026）
已經(jīng)有0人回復(fù)
-
物理學(xué)I論文潤色/翻譯怎么收費(fèi)?
已經(jīng)有103人回復(fù)
-
哈爾濱理工大學(xué)物理系招收物理學(xué)考研調(diào)劑
已經(jīng)有5人回復(fù)
-
0702一志愿吉大B區(qū)求調(diào)劑
已經(jīng)有5人回復(fù)
-
求調(diào)劑
已經(jīng)有0人回復(fù)
-
0702一志愿吉大B區(qū)求調(diào)劑有論文
已經(jīng)有0人回復(fù)
-
請問還有沒有用Latex寫文章的小伙伴們？
已經(jīng)有0人回復(fù)
-
光學(xué)工程學(xué)碩調(diào)劑信息
已經(jīng)有30人回復(fù)
-
歡迎加入課題組
已經(jīng)有0人回復(fù)
-
散金幣，求好運(yùn)，祝面上順利！
已經(jīng)有61人回復(fù)
dingd
鐵桿木蟲 (職業(yè)作家)
- 計算強(qiáng)帖: 4
- 應(yīng)助: 1641 (講師)
- 金幣: 15037.3
- 散金: 101
- 紅花: 234
- 帖子: 3410
- 在線: 1223.7小時
- 蟲號: 291104
- 注冊: 2006-10-28
2樓2017-10-11 09:19:32
manmanbobo
新蟲 (初入文壇)
- 應(yīng)助: 0 (幼兒園)
- 金幣: 25
- 紅花: 1
- 帖子: 37
- 在線: 12.5小時
- 蟲號: 2965461
- 注冊: 2014-02-09
- 專業(yè): 地理信息系統(tǒng)
3樓2017-10-16 10:31:21
相關(guān)版塊跳轉(zhuǎn)
我要訂閱樓主 manmanbobo 的主題更新
| 3 | 1/1 | 返回列表 |
| 最具人氣熱帖推薦 [查看全部] | 作者 | 回/看 | 最后發(fā)表 | |
|---|---|---|---|---|
|
[考研] 321求調(diào)劑 +8 | 何潤采123 2026-03-18 | 10/500 |
|
|---|---|---|---|---|
|
[考博] 東華理工大學(xué)化材專業(yè)26屆碩士博士申請 +8 | zlingli 2026-03-13 | 8/400 |
|
|
|
[考研] 求調(diào)劑，一志愿:南京航空航天大學(xué)大學(xué) ，080500材料科學(xué)與工程學(xué)碩，總分289分 +3 | @taotao 2026-03-19 | 3/150 |
|
|
|
[教師之家] 焦慮 +9 | 水冰月月野兔 2026-03-13 | 13/650 |
|
|
|
[考研] 274求調(diào)劑 +6 | S.H1 2026-03-18 | 6/300 |
|
|
|
[考研] 354求調(diào)劑 +4 | Tyoumou 2026-03-18 | 7/350 |
|
|
|
[考研] 311求調(diào)劑 +6 | 26研0 2026-03-15 | 6/300 |
|
|
|
[考研] 070300化學(xué)319求調(diào)劑 +6 | 錦鯉0909 2026-03-17 | 6/300 |
|
|
|
[考研] 工科材料085601 279求調(diào)劑 +6 | 困于星晨 2026-03-17 | 6/300 |
|
|
|
[考研] 301求調(diào)劑 +9 | yy要上岸呀 2026-03-17 | 9/450 |
|
|
|
[考研] 296求調(diào)劑 +5 | 大口吃飯 身體健 2026-03-13 | 5/250 |
|
|
|
[考研] 考研化學(xué)學(xué)碩調(diào)劑，一志愿985 +4 | 張vvvv 2026-03-15 | 6/300 |
|
|
|
[考研] 290求調(diào)劑 +3 | p asserby. 2026-03-15 | 4/200 |
|
|
|
[考研] 302求調(diào)劑 +4 | 小賈同學(xué)123 2026-03-15 | 8/400 |
|
|
|
[考研] 274求調(diào)劑 +5 | 時間點(diǎn) 2026-03-13 | 5/250 |
|
|
|
[考研] 東南大學(xué)364求調(diào)劑 +5 | JasonYuiui 2026-03-15 | 5/250 |
|
|
|
[考研] 一志愿211 0703方向310分求調(diào)劑 +3 | 努力奮斗112 2026-03-15 | 3/150 |
|
|
|
[考研] 070303 總分349求調(diào)劑 +3 | LJY9966 2026-03-15 | 5/250 |
|
|
|
[考研] 材料與化工 323 英一+數(shù)二+物化，一志愿：哈工大 本人本科雙一流 +4 | 自由的_飛翔 2026-03-13 | 5/250 |
|
|
|
[考研] 308求調(diào)劑 +3 | 是Lupa啊 2026-03-12 | 3/150 |
|