我将一个非线性模型与一个数据集进行了拟合.
但是，我需要执行该模型的半正常图(使用hnp包中的hnp()函数).
考虑到stats包中的nls()函数，我的模型是适合的，请参见:
nonlinear_func <- function(x, beta_0, beta_1) {
return (beta_0 * (1 - exp(beta_1 * x)) + 1.30)
}
nonlinear_model <- nls(y ~ nonlinear_func(x, beta_0, beta_1),
start = list(beta_0 = 16, beta_1 = -0.2),
data.set)
数据集如下:
data.set <- structure(list(x = c(13.05, 6.05, 13.21, 9.55, 18.14, 9.55, 14.48,
15.28, 9.87, 15.92, 12.41, 12.41, 12.57, 15.12, 10.66, 16.87,
12.57, 15.92, 9.71, 15.92, 17.35, 6.37, 11.94, 11.14, 8.91, 13.05,
17.67, 10.66, 17.19, 7, 10.82, 11.62, 16.71, 18.3, 11.78, 12.89,
10.82, 9.23, 14.32, 7.64, 5.09, 15.44, 10.35, 8.91, 14.32, 13.21,
8.91, 15.6, 14.16, 15.28, 12.57), y = c(16.8, 11.3, 16.9, 11.8,
17.5, 15.8, 20, 19.1, 15.8, 18.8, 18.6, 18.4, 18.4, 18.7, 15.2,
18.3, 15.7, 17.5, 14.8, 16.7, 19.8, 10.3, 18.6, 14.4, 14.3, 17.8,
21, 18.8, 18.9, 13.7, 17.6, 17.5, 19.6, 18.8, 15.3, 17, 15.9,
13.3, 17.3, 13.6, 9.3, 17.7, 14.2, 14.9, 18.4, 18.2, 14.3, 19.7,
18.6, 18.1, 15.5)), row.names = c(17L, 18L, 19L, 20L, 21L, 22L,
23L, 24L, 25L, 26L, 27L, 28L, 29L, 30L, 31L, 32L, 33L, 34L, 35L,
36L, 37L, 38L, 39L, 40L, 41L, 42L, 43L, 44L, 45L, 46L, 47L, 48L,
49L, 50L, 51L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 11L,
12L, 13L, 14L, 15L, 16L), class = "data.frame")
好的，当从hnp包(https://cran.r-project.org/web/packages/hnp/hnp.pdf)第10页、第14页和第15页读取信息时，您可以看到没有与nls()函数相关联的实现，但是可以通过几个步骤来构建图，例如，参见使用gamlss()中的调整模型的示例(hnp包中没有实现的另一个函数):
## Example no. 2: Implementing gamma model using package gamlss
# load package
library(gamlss)
# model fitting
y <- rGA(30, mu=rep(c(.5, 1.5, 5), each=10), sigma=.5)
tr <- gl(3, 10)
fit2 <- gamlss(y ~ tr, family=GA)
# diagfun
d.fun <- function(obj) resid(obj) # this is the default if no
# diagfun is provided
# simfun
s.fun <- function(n, obj) {
mu <- obj$mu.fv
sig <- obj$sigma.fv
rGA(n, mu=mu, sigma=sig)
}
# fitfun
my.data <- data.frame(y, tr)
f.fun <- function(y.) gamlss(y. ~ tr, family=GA, data=my.data)
# hnp call
hnp(fit2, newclass=TRUE, diagfun=d.fun, simfun=s.fun,
fitfun=f.fun, data=data.frame(y, tr))
考虑到这些信息，我try 对我的非线性模型执行同样的操作，请参见:
d.fun <- function(obj) resid(obj)
s.fun <- function(n, obj) {
}
f.fun <- function(data) {
nls(y ~ nonlinear_func(x, beta_0, beta_1),
start = list(beta_0 = 16, beta_1 = -0.2),
data = data.set)
}
library(hnp)
hnp(nonlinear_model, newclass = TRUE,
diagfun = d.fun, simfun = s.fun,
fitfun = f.fun, data = data.set)
但是，请注意发生了一些错误，因为在创建图形时没有考虑模拟包络，请参见:
在car包中有一个名为qqPlot的函数，它执行下面的图(这不是半正常的曲线图)，但是我需要独占使用hnp包.
library(car)
qqPlot(residuals(nonlinear_model), pch = 19, col = "green", cex = 1.2,
xlab = "Norm Quantiles",
ylab = "Residuals")
