Начальные значения для 4-х параметров NLS - функция Чепмена-Ричардса - PullRequest
Новая
       19

Начальные значения для 4-х параметров NLS - функция Чепмена-Ричардса

0 голосов
/ 20 ноября 2018

* Примечание - я прочитал несколько постов о том, как найти начальные значения для NLS - однако я не нашел ни одного с уравнением этой формы (т.е. 4 параметра, показатель степени возведен в степень)

Я изо всех сил пытаюсь найти подходящие начальные значения для уравнения Чепмена-Ричардса, которое обычно используется в лесном хозяйстве для моделирования роста деревьев. 

y(t) = α * (1 - β * exp(-k * t)^{1/(1-m)})

Обычно я пытаюсь найти начальные значения путем построения линии с наборомпараметров, а затем настраивая их так, чтобы они более точно подходили к данным (Рисунок 1).После этого я бы использовал параметры в функции: 

initial.test <- chapmanRichards(seq(0:15),42,0.95,0.28, 0.67)
plot(age,topHeight,type="p",xlab="year since planting",ylab="Dom height (m)", xlim = c(0,20), ylim = c(0, 50))
lines(seq(0:15),initial.test,col="red")

enter image description here

nls(topHeight ~ chapmanRichards(age,a,b,k,m),start=list(a=42,b=0.95,k=0.28,m=0.67))

В этом случае программа может соответствовать кривойс предоставленными начальными значениями.Проблема, однако, заключается в том, что когда данные немного зашумлены, и после 2 часов работы с начальными значениями теста я все еще не могу найти достаточно хороших начальных значений (на рисунке 2 показаны несколько попыток для другого набора данных.

enter image description here

Кто-нибудь может посоветовать, как можно найти подходящие начальные значения? Я думал о создании матрицы, которая в основном запускает последовательность для каждого изпараметры и цикл nls с этими начальными значениями, но не уверен, как будет выглядеть код. Любой другой совет будет с благодарностью!

PS - это будет что-то более подходящее для Excel - решатель?

1 Ответ

0 голосов
/ 22 ноября 2018

Как отметил @Roland в комментариях, параметры в уравнении, показанном в вопросе, не могут быть идентифицированы, поэтому при условии, что уравнение такое, как он показал: 

y = a * (1 - b * exp(-k * t))^{1/(1-m)}

взять журнал обеих сторон: 

log(y) ~ log(a) + (1/(1-m)) * log(1 - b * exp(-k*t))

и пусть log (a) = A, 1 / (1-m) = M и b = exp (k * B) дает: 

log(y) ~ A + M * log(1 - exp(k*(B-t))

Так как B является смещениеми k - это масштабирование, мы можем оценить их как B = mean (t) и k = 1 / sd (t).Используя algorithm = "plinear", мы можем избежать начальных значений для линейных параметров (A и M), при условии, что мы указываем правую часть в виде матрицы таким образом, чтобы A, умноженное на первый столбец, плюс M, умноженное на второй столбец, давало прогнозируемое значение.Таким образом, мы имеем: 

st <- list(B = mean(t), k = 1/sd(t))
fm0 <- nls(log(y) ~ cbind(1, log(1 - exp(k*(B - t)))), start = st,
  algorithm = "plinear")

и затем обратно преобразовываем полученные коэффициенты, чтобы получить начальные значения для запуска окончательного nls.

Также обратите внимание, что nls2 в пакете nls2можно оценить модель по сетке или по случайному набору точек, чтобы получить начальные значения.

...