I think that these aren't good initial values. If you do a plot of data and
add a curve, the curve don't approximate the data. Frequently I use
interactive procedures to get good initial values. Using playwith() you can
handle sliders to adjust values and use in nls(), look the following
r <-
c(1.16,1.143,1.109,1.093,1.079,1.066,1.053,1.040,1.027,1.015,1.004,0.994,0.985,0.977)
D <-
c(0.1806551,0.2703113,0.3757225,0.5271811,0.8665835,1.0812568,1.0612762,1.0726612,
1.167927,1.191092,1.1336938,1.1215107,0.9619603,0.8315467)
names <- c(La,Ce,Pr,Nd,Sm,Eu,Gd,Tb,Dy,Ho,Er,Tm,Yb,Lu)
plot(log(D)~r)
Do <- 0.8; En <- 390; ro <- 1.03
curve(log10(Do)*exp(((-4*pi*En*Na)*((ro/2)*(x-ro)^2+(1/3)*(x-ro)^3))/(R*T)),
add=TRUE)
Na <- 6.0221415*10^23 # Avrogadro's number
T <- 1010 # Temp in K
R <- 8.3144 # Gas constant [J mol^-1 K^-1]
## Brice-Model
bricemod <-
nls(log10(D)~log10(Do)*exp(((-4*pi*En*Na)*((ro/2)*(r-ro)^2+(1/3)*(r-ro)^3))/(R*T)),
start=list(Do=0.8, En=390, ro=1.03), trace=TRUE)
require(playwith)
start <- list()
playwith(
{
plot(log(D)~r)
curve(log10(Do)*exp(((-4*pi*En*Na)*((ro/2)*(x-ro)^2+(1/3)*(x-ro)^3))/(R*T)),
add=TRUE)
start <<- list(Do=Do, En=En, ro=ro)
},
parameters=list(
Do=seq(0.7,0.8,by=0.001),
En=seq(350,450,by=1),
ro=seq(0.5,1.5,by=0.1)))
start
bricemod <-
nls(log10(D)~log10(Do)*exp(((-4*pi*En*Na)*((ro/2)*(r-ro)^2+(1/3)*(r-ro)^3))/(R*T)),
start=start, trace=TRUE)
Bests.
Walmes.
==========================================================================
Walmes Marques Zeviani
LEG (Laboratório de Estatística e Geoinformação, 25.450418 S, 49.231759 W)
Departamento de Estatística - Universidade Federal do Paraná
fone: (+55) 41 3361 3573
VoIP: (3361 3600) 1053 1173
e-mail: wal...@ufpr.br
twitter: @walmeszeviani
homepage: http://www.leg.ufpr.br/~walmes
linux user number: 531218
==========================================================================
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
