topic Re: Issues with Convergence in Statistical Procedures

Issues with Convergence

Nebulus — Fri, 19 Sep 2014 00:34:59 GMT

So I've been trying so solve this problem for a few weeks and am at my whit's end at this point. I've tried all the troubleshooting I can to no avail. Whenever I try to run this model (rate= a*temp*(Temp-Tmin)*((Tmax-temp)**(1/2)) it always fails to converge. The issue seems to be with (Tmax-temp)**(1/2). Ive tried Sqrt(Tmax-temp) and (Tmax-temp)**(0.5) as well.

Here is my code with my data if someone wants to try and see for themselves what happens.

*options ls=78 ps=55 nocenter nodate nonumber;

title "1st Instar Female Development";

data dev1;

input id days rate temp;

if temp = 34.9 then delete;

datalines;

99 0 0 34.9

100 0 0 34.9

101 0 0 29.54

102 0 0 29.54

113 0 0 29.54

114 0 0 29.54

115 0 0 29.54

116 0 0 29.54

117 0 0 29.54

118 0 0 29.54

119 0 0 29.54

120 0 0 29.54

121 0 0 29.54

122 0 0 29.54

123 0 0 29.54

124 0 0 29.54

125 0 0 29.54

126 0 0 29.54

127 0 0 29.54

128 0 0 29.54

129 0 0 29.54

130 0 0 29.54

131 0 0 29.54

132 0 0 29.54

133 0 0 29.54

134 0 0 29.54

135 0 0 29.54

136 0 0 29.54

137 0 0 29.54

138 0 0 29.54

139 0 0 29.54

140 0 0 29.54

141 0 0 29.54

142 0 0 29.54

143 0 0 29.54

144 0 0 29.54

145 0 0 29.54

146 0 0 29.54

147 0 0 29.54

148 0 0 29.54

149 0 0 29.54

150 0 0 29.54

151 0 0 29.54

152 0 0 29.54

153 0 0 29.54

154 0 0 29.54

155 0 0 29.54

156 0 0 29.54

157 0 0 29.54

158 0 0 29.54

159 0 0 29.54

160 0 0 29.54

161 0 0 29.54

162 0 0 29.54

163 0 0 29.54

164 0 0 29.54

165 0 0 29.54

166 0 0 29.54

167 0 0 29.54

168 0 0 29.54

169 0 0 29.54

170 0 0 29.54

171 0 0 29.54

172 0 0 29.54

173 0 0 29.54

174 0 0 29.54

175 0 0 29.54

176 0 0 29.54

177 0 0 29.54

178 0 0 29.54

179 0 0 29.54

180 0 0 29.54

181 0 0 29.54

182 0 0 29.54

183 0 0 29.54

184 0 0 29.54

185 0 0 29.54

186 0 0 29.54

187 0 0 29.54

188 0 0 29.54

189 0 0 29.54

190 0 0 29.54

191 0 0 29.54

192 0 0 29.54

193 0 0 29.54

194 0 0 29.54

195 0 0 29.54

196 0 0 29.54

197 0 0 29.54

198 0 0 29.54

199 0 0 29.54

200 0 0 29.54

201 9.45 0.105820106 27.91

202 3.95 0.253164557 27.91

203 3.95 0.253164557 27.91

204 4.95 0.202020202 27.91

205 3.95 0.253164557 27.91

206 3.95 0.253164557 27.91

207 3.95 0.253164557 27.91

208 4.95 0.202020202 27.91

209 4.95 0.202020202 27.91

210 9.75 0.102564103 27.91

211 4.75 0.210526316 27.91

212 4.75 0.210526316 27.91

213 5.75 0.173913043 27.91

214 4.75 0.210526316 27.91

215 3.75 0.266666667 27.91

216 4.75 0.210526316 27.91

217 5.4 0.185185185 27.91

218 5.4 0.185185185 27.91

219 3.4 0.294117647 27.91

220 5.4 0.185185185 27.91

221 4.4 0.227272727 27.91

222 4.39 0.227790433 27.91

223 4.06 0.246305419 27.91

224 4.06 0.246305419 27.91

225 5.88 0.170068027 27.91

226 6.88 0.145348837 27.91

227 4.88 0.204918033 27.91

228 5.88 0.170068027 27.91

229 2.88 0.347222222 27.91

230 4.88 0.204918033 27.91

231 5.88 0.170068027 27.91

232 4.75 0.210526316 27.91

233 4.75 0.210526316 27.91

234 4.75 0.210526316 27.91

235 5.75 0.173913043 27.91

236 5.75 0.173913043 27.91

237 4.75 0.210526316 27.91

238 2.75 0.363636364 27.91

239 4.75 0.210526316 27.91

240 2.75 0.363636364 27.91

241 3.75 0.266666667 27.91

242 3.75 0.266666667 27.91

243 3.9375 0.253968254 25.3

244 2.9375 0.340425532 25.3

245 3.9375 0.253968254 25.3

246 3.9375 0.253968254 25.3

247 2.9375 0.340425532 25.3

248 3.9375 0.253968254 25.3

249 2.9375 0.340425532 25.3

250 2.9375 0.340425532 25.3

251 3.9375 0.253968254 25.3

252 2.9375 0.340425532 25.3

253 3.9375 0.253968254 25.3

254 3.9375 0.253968254 25.3

255 4.9375 0.202531646 25.3

256 4.9375 0.202531646 25.3

257 4.9375 0.202531646 25.3

258 4.9375 0.202531646 25.3

259 4.9375 0.202531646 25.3

260 3.9375 0.253968254 25.3

261 3.9375 0.253968254 25.3

262 6.9375 0.144144144 25.3

263 4 0.25 25.3

264 4 0.25 25.3

265 4 0.25 25.3

266 3 0.333333333 25.3

267 5 0.2 25.3

268 4 0.25 25.3

269 5 0.2 25.3

270 5 0.2 25.3

271 5 0.2 25.3

272 5 0.2 25.3

273 5 0.2 25.3

274 5 0.2 25.3

275 4.42 0.226244344 25.3

276 3.42 0.292397661 25.3

277 4.42 0.226244344 25.3

278 2.42 0.41322314 25.3

279 3.42 0.292397661 25.3

280 6.42 0.15576324 25.3

281 5.95 0.168067227 25.3

282 6.95 0.143884892 25.3

283 3.95 0.253164557 25.3

284 3.95 0.253164557 25.3

285 2.95 0.338983051 25.3

286 3.95 0.253164557 25.3

287 3.95 0.253164557 25.3

288 3.95 0.253164557 25.3

289 2.95 0.338983051 25.3

290 3.95 0.253164557 25.3

291 4.95 0.202020202 25.3

292 4.95 0.202020202 25.3

293 3.95 0.253164557 25.3

294 3.95 0.253164557 25.3

295 8.69 0.115074799 20.24

296 11.69 0.085543199 20.24

297 6.69 0.149476831 20.24

298 11.69 0.085543199 20.24

299 7.69 0.130039012 20.24

300 6.53 0.153139357 20.24

301 6.53 0.153139357 20.24

302 9.53 0.104931794 20.24

303 7.53 0.132802125 20.24

304 7.53 0.132802125 20.24

305 7.53 0.132802125 20.24

306 6.53 0.153139357 20.24

307 7.53 0.132802125 20.24

308 6.53 0.153139357 20.24

309 5.83 0.171526587 20.24

310 6.83 0.146412884 20.24

311 5.83 0.171526587 20.24

312 4.83 0.207039337 20.24

313 5.83 0.171526587 20.24

314 9.83 0.1017294 20.24

315 8.83 0.113250283 20.24

316 7.83 0.127713921 20.24

317 6.83 0.146412884 20.24

318 8.83 0.113250283 20.24

319 6.83 0.146412884 20.24

320 5.83 0.171526587 20.24

321 8.83 0.113250283 20.24

322 11.63 0.085984523 20.24

323 10.63 0.094073377 20.24

324 6.63 0.150829563 20.24

325 7.63 0.131061599 20.24

326 7.63 0.131061599 20.24

327 6.63 0.150829563 20.24

328 6.63 0.150829563 20.24

329 12.63 0.079176564 20.24

330 6.63 0.150829563 20.24

331 6.63 0.150829563 20.24

332 15.95 0.062695925 14.8

333 12.95 0.077220077 14.8

334 13.95 0.071684588 14.8

335 19.95 0.050125313 14.8

336 16.95 0.05899705 14.8

337 12.95 0.077220077 14.8

338 17.95 0.055710306 14.8

339 15.95 0.062695925 14.8

340 13.95 0.071684588 14.8

341 24.95 0.04008016 14.8

342 12.95 0.077220077 14.8

343 17.58 0.056882821 14.8

344 23.58 0.042408821 14.8

345 21.58 0.046339203 14.8

346 24.58 0.040683483 14.8

347 21.58 0.046339203 14.8

348 18.58 0.053821313 14.8

349 19.58 0.051072523 14.8

350 15.58 0.064184852 14.8

351 17.67 0.056593096 14.8

352 20.67 0.048379294 14.8

353 14.67 0.068166326 14.8

354 18.67 0.053561864 14.8

355 19.67 0.050838841 14.8

356 15.67 0.063816209 14.8

357 18.67 0.053561864 14.8

358 18 0.055555556 14.8

359 18 0.055555556 14.8

360 17 0.058823529 14.8

361 20 0.05 14.8

362 14.64 0.068306011 14.8

363 15.64 0.063938619 14.8

364 15.64 0.063938619 14.8

365 14.64 0.068306011 14.8

366 17.64 0.056689342 14.8

367 16.64 0.060096154 14.8

368 14.64 0.068306011 14.8

369 16.64 0.060096154 14.8

370 16.53 0.060496068 14.8

371 17.53 0.057045066 14.8

372 19.53 0.051203277 14.8

373 22.53 0.044385264 14.8

374 27 0.037037037 14.8

375 16 0.0625 14.8

376 14 0.071428571 14.8

377 17 0.058823529 14.8

378 19 0.052631579 14.8

379 20 0.05 14.8

380 10.42 0.09596929 14.8

381 16.42 0.06090134 14.8

382 20.42 0.048971596 14.8

383 15.42 0.064850843 14.8

384 15.42 0.064850843 14.8

385 15.42 0.064850843 14.8

386 17.42 0.057405281 14.8

387 14.42 0.069348128 14.8

388 15.42 0.064850843 14.8

389 22.42 0.044603033 14.8

390 15.74 0.063532402 14.8

391 0 0 10.17

392 0 0 10.17

393 0 0 10.17

394 0 0 10.17

395 0 0 10.17

396 0 0 10.17

397 0 0 10.17

398 0 0 10.17

399 0 0 10.17

400 0 0 10.17

401 0 0 10.17

402 0 0 10.17

403 0 0 10.17

404 0 0 10.17

405 0 0 10.17

406 0 0 10.17

407 0 0 10.17

408 0 0 10.17

409 0 0 10.17

410 0 0 10.17

411 0 0 10.17

412 0 0 10.17

413 0 0 10.17

414 0 0 10.17

415 0 0 10.17

416 0 0 10.17

417 0 0 10.17

418 0 0 10.17

419 0 0 10.17

420 0 0 10.17

421 0 0 10.17

422 0 0 10.17

423 0 0 10.17

424 0 0 10.17

425 0 0 10.17

426 0 0 10.17

427 0 0 10.17

428 0 0 10.17

429 0 0 10.17

430 0 0 10.17

431 0 0 10.17

432 0 0 10.17

433 0 0 10.17

434 0 0 10.17

435 0 0 10.17

436 0 0 10.17

437 0 0 10.17

438 0 0 10.17

439 0 0 10.17

440 0 0 10.17

441 0 0 10.17

442 0 0 10.17

443 0 0 10.17

444 0 0 10.17

445 0 0 10.17

446 0 0 10.17

447 0 0 10.17

448 0 0 10.17

449 0 0 10.17

450 0 0 10.17

451 0 0 10.17

452 0 0 10.17

453 0 0 10.17

454 0 0 10.17

455 0 0 10.17

456 0 0 10.17

457 0 0 10.17

458 0 0 10.17

459 0 0 10.17

460 0 0 10.17

461 0 0 10.17

462 0 0 10.17

463 0 0 10.17

464 0 0 10.17

465 0 0 10.17

466 0 0 10.17

467 0 0 10.17

468 0 0 10.17

469 0 0 10.17

470 0 0 10.17

471 0 0 10.17

472 0 0 10.17

473 0 0 10.17

474 0 0 10.17

475 0 0 10.17

476 0 0 10.17

477 0 0 10.17

478 0 0 10.17

479 0 0 10.17

480 0 0 10.17

481 0 0 10.17

482 0 0 10.17

483 0 0 10.17

484 0 0 10.17

485 0 0 10.17

486 0 0 10.17

487 0 0 10.17

488 0 0 10.17

489 0 0 10.17

490 0 0 10.17

;

data DevoFiller;

do temp = 0 to 35 by 0.05;

predict=1;

output;

end;

run;

data dev;

set dev1 DevoFiller;

run;

proc sort data=dev;

by temp;

quit;

proc print data=dev;

quit;

*Title1 'Setup Marquardt, a, Tmin, and Tmax, No Briere Derivatives,';

proc nlin data=dev

Best=10

method=Marquardt /*Need to specify the algorithm that will be used to estimate the four parameters.

Apparently you use the Levenberg-Marquardt method when parameters are highly correlated*/

converge=0.00001 /*Specifying the convergence criterion listed by Lactin et al.*/

; /*Specifying the maximum number of iterations to estimate parameters*/

parms /*Specifying range of values for parameters to conduct grid search*/

a= 0.0001 to 0.001 by 0.001

Tmin= 10 to 15 by 0.1

Tmax= 28 to 35 by 0.1;

model rate = a*temp*(temp-Tmin)*((Tmax-temp)*(1/2));

/*Model outlined by Briere et al. 1999*/

output out=Devonlinout predicted=predDevo L95M=Cl95meanDevo U95M=CU95meanDevo L95=CL95indDevo U95=CU95indDevo;

quit;

proc print data=Devonlinout;

quit;

proc means data=dev CSS;

var rate;

quit;

Any advice would be greatly useful!

Re: Issues with Convergence

PGStats — Fri, 19 Sep 2014 02:37:22 GMT

You must prevent parameters tmin and tmax from wandering into the data range for this to converge. I got a fit in 4 iterations using default methods this way:

proc sql;

select

min(temp)- 0.001,

max(temp) + 0.001,

max(rate) / (mean(temp)*(mean(temp)-min(temp))*sqrt(max(temp)-mean(temp)))

into :tmin, :tmax, :a from dev1;

quit;

data dev2;

set dev1 end=done;

if done then do;

call missing(rate);

predict + 1;

do temp = &tmin to &tmax by 0.05;

output;

end;

else output;

run;

proc sort data=dev2; by temp; run;

proc nlin data=dev2;

parms tmin=&tmin, tmax=&tmax, a=&a;

bounds tmin <= &tmin, tmax >= &tmax;

model rate = a * temp * (temp-tmin) * sqrt(tmax-temp);

output out=dev3 predicted=predRate;

run;

data devgraph;

set dev3;

if not missing(rate) then call missing(predRate);

run;

ods listing gpath="&sasforum.\graphs";

ods graphics / imagename="instar" reset=index;

proc sgplot data=devGraph;

scatter x=temp y=rate;

series x=temp y=predRate;

run;

Re: Issues with Convergence

Nebulus — Fri, 19 Sep 2014 04:15:11 GMT

Thanks a lot. So it would seem that it will always predict Tmax and Tmin at exactly my observed Tmin and Tmax? This seems to be a bit different that what other literature is getting for predicted vs obs. tmin/tmax though. Looking at the predicted rate corresponding to the predicted Tmax it is still a positive value and not near 0 yet.

I haven't worked with bounds before.. Is their table values a correction to their corresponding prediction?

Parameter Estimate ApproxStd Error Appro. 95% CL | Label|

tmin 10.1690 0 10.1690 /10.1690

tmax 29.5410 0 29.5410 /29.5410

a 0.000297 4.372E-6 0.000288 /0.000305

Bound0 0.0268 0.00579 0.0155 /0.0382 | tmin <= 10.169 |

Bound1 1.3629 0.9638 -0.5250 /3.2508 |29.541 <= tmax |

So for example Tmin= 10.169+ 0.0268 and Tmax= 29.54+1.3629 with the SE listed?

The paper this model came from states that Temp <= Tlow and T >= Tmax. That is essentially what you did correct?

Re: Issues with Convergence

PGStats — Fri, 19 Sep 2014 21:59:39 GMT

Actually, looking at your data, it would be reasonable to assume that tmin is somewhere between 11 and 14.5 and tmax is somewhere between 28 and 29.5 . But to allow these parameter values, we must modify the model so that predicted rates are zero outside the range (tmin-tmax).

Also, you may want to explore other curvatures for the upward limb of the curve. For this, I added the exponent parameter. It yields decent models for values in the range 0.6 - 1.45.

%let exponent=1.0;

proc sql;

select

min(temp)- 0.001,

max(temp) + 0.001,

max(rate) / (mean(temp)*(mean(temp)-min(temp))**(&exponent)*sqrt(max(temp)-mean(temp)))

into :tmin0, :tmax0, :a0

from dev1;

select

min(temp)- 0.001,

max(temp) + 0.001,

max(rate) / (mean(temp)*(mean(temp)-min(temp))**(&exponent)*sqrt(max(temp)-mean(temp)))

into :tmin, :tmax, :a

from dev1

where rate > 0.001;

quit;

data dev2;

set dev1 end=done;

if done then do;

call missing(rate);

predict + 1;

do temp = &tmin0 to &tmax0 by 0.05, &tmax0;

output;

end;

else output;

run;

proc sort data=dev2; by temp; run;

proc nlin data=dev2 convergeparm=1e-7;

parms tmin=&tmin0, tmax=&tmax0, a=&a;

bounds &tmin0 <= tmin < &tmin, &tmax < tmax <= &tmax0;

if temp <= tmin then r = 0;

else if temp >= tmax then r = 0;

else r = a * temp * (temp-tmin)**(&exponent) * sqrt(tmax-temp);

model rate = r;

output out=dev3 predicted=predRate;

run;

data devgraph;

set dev3;

if not missing(rate) then call missing(predRate);

run;

ods listing gpath="&sasforum.\graphs";

ods graphics / imagename="instar-exponent &exponent." reset=index;

title2 "Exponent=&exponent.";

proc sgplot data=devGraph;

scatter x=temp y=rate;

series x=temp y=predRate;

run;

Re: Issues with Convergence

Nebulus — Sat, 20 Sep 2014 02:25:58 GMT

Thanks so much. You actually hit the nail on the head.

The model I originally posted is the Briere-1 Model: aT(T-Tmin)(Tmax-Tmin)**(1/2)

He also proposed a more generalized model Briere-2: aT(T-Tmin)(Tmax-Tmin)**(1/m), where m is a 4th parameter to estimate.

Using your code and adding m I predicted m ~ 2.75. By adding this term the RSS went from 0.500 to 0.4096. So yes by modifying the exponent you tend to get a better fit of the data.