R 软件 计算也很方便,例子如下,供参考
http://127.0.0.1:20985/library/nlsur/html/nlsur.html,这是帮助文件
##
https://janmarvin.r-universe.dev/nlsur
## install.packages('nlsur', repos = c('https://janmarvin.r-universe.dev', 'https://cloud.r-project.org'))
library(nlsur)
data(costs)
dat <- costs
# apply a patch to create Greene Ed. 7 Data
dat$Sm[dat$Year == 1958] <- 0.61886
dat$Pe[dat$Year == 1950] <- 1.12442
dat$Pm[dat$Year == 1949] <- 1.06625
# model
model <- list(
Sk ~ bk + dkk * log(Pk/Pm) + dkl * log(Pl/Pm) + dke * log(Pe/Pm),
Sl ~ bl + dkl * log(Pk/Pm) + dll * log(Pl/Pm) + dle * log(Pe/Pm),
Se ~ be + dke * log(Pk/Pm) + dle * log(Pl/Pm) + dee * log(Pe/Pm)
)
erg <- nlsur(eqns = model, data = dat, type = "FGNLS")
erg
## Not run:
# Greene Example 10.3
library(nlsur)
url <- "http://www.stern.nyu.edu/~wgreene/Text/Edition7/TableF10-2.txt"
dd <- read.table(url, header = T)
names(dd) <-
c("Year", "Cost", "Sk", "Sl", "Se", "Sm", "Pk", "Pl", "Pe", "Pm")
eqns <-
list( Sk ~ bk + dkk * log(Pk/Pm) + dkl * log(Pl/Pm) + dke * log(Pe/Pm),
Sl ~ bl + dkl * log(Pk/Pm) + dll * log(Pl/Pm) + dle * log(Pe/Pm),
Se ~ be + dke * log(Pk/Pm) + dle * log(Pl/Pm) + dee * log(Pe/Pm))
strtvls <- c(be = 0, bk = 0, bl = 0,
dkk = 0, dkl = 0, dke = 0,
dll = 0, dle = 0, dee = 0)
erg <- ?nlsur(eqns = eqns, data = dd, startvalues = strtvls, type ="FGNLS",
trace = TRUE, eps = 1e-10)
erg