- html - 出于某种原因,IE8 对我的 Sass 文件中继承的 html5 CSS 不友好?
- JMeter 在响应断言中使用 span 标签的问题
- html - 在 :hover and :active? 上具有不同效果的 CSS 动画
- html - 相对于居中的 html 内容固定的 CSS 重复背景?
我使用smooth.spline
来估计我的数据的三次样条。但是当我使用方程计算 90% 逐点置信区间时,结果似乎有点偏差。有人可以告诉我我是否做错了吗?我只是想知道是否有一个函数可以自动计算与 smooth.spline
函数关联的逐点间隔带。
boneMaleSmooth = smooth.spline( bone[males,"age"], bone[males,"spnbmd"], cv=FALSE)
error90_male = qnorm(.95)*sd(boneMaleSmooth$x)/sqrt(length(boneMaleSmooth$x))
plot(boneMaleSmooth, ylim=c(-0.5,0.5), col="blue", lwd=3, type="l", xlab="Age",
ylab="Relative Change in Spinal BMD")
points(bone[males,c(2,4)], col="blue", pch=20)
lines(boneMaleSmooth$x,boneMaleSmooth$y+error90_male, col="purple",lty=3,lwd=3)
lines(boneMaleSmooth$x,boneMaleSmooth$y-error90_male, col="purple",lty=3,lwd=3)
因为我不确定我是否正确执行,所以我使用了 mgcv
包中的 gam()
函数。
它立即给出了一个置信区间,但我不确定它是 90% 还是 95% CI 还是其他。如果有人能解释一下就太好了。
males=gam(bone[males,c(2,4)]$spnbmd ~s(bone[males,c(2,4)]$age), method = "GCV.Cp")
plot(males,xlab="Age",ylab="Relative Change in Spinal BMD")
最佳答案
我不确定 smooth.spline
的置信区间是否像 lowess
那样具有“良好”的置信区间。但我从 CMU Data Analysis course 找到了一个代码示例建立贝叶斯引导置信区间。
以下是使用的函数和示例。主要函数是 spline.cis,其中第一个参数是数据框,其中第一列是 x 值,第二列是 y 值> 值(value)观。另一个重要参数是 B,它指示要执行的引导复制数量。 (有关完整详细信息,请参阅上面链接的 PDF。)
# Helper functions
resampler <- function(data) {
n <- nrow(data)
resample.rows <- sample(1:n,size=n,replace=TRUE)
return(data[resample.rows,])
}
spline.estimator <- function(data,m=300) {
fit <- smooth.spline(x=data[,1],y=data[,2],cv=TRUE)
eval.grid <- seq(from=min(data[,1]),to=max(data[,1]),length.out=m)
return(predict(fit,x=eval.grid)$y) # We only want the predicted values
}
spline.cis <- function(data,B,alpha=0.05,m=300) {
spline.main <- spline.estimator(data,m=m)
spline.boots <- replicate(B,spline.estimator(resampler(data),m=m))
cis.lower <- 2*spline.main - apply(spline.boots,1,quantile,probs=1-alpha/2)
cis.upper <- 2*spline.main - apply(spline.boots,1,quantile,probs=alpha/2)
return(list(main.curve=spline.main,lower.ci=cis.lower,upper.ci=cis.upper,
x=seq(from=min(data[,1]),to=max(data[,1]),length.out=m)))
}
#sample data
data<-data.frame(x=rnorm(100), y=rnorm(100))
#run and plot
sp.cis <- spline.cis(data, B=1000,alpha=0.05)
plot(data[,1],data[,2])
lines(x=sp.cis$x,y=sp.cis$main.curve)
lines(x=sp.cis$x,y=sp.cis$lower.ci, lty=2)
lines(x=sp.cis$x,y=sp.cis$upper.ci, lty=2)
这给出了类似的东西
实际上,看起来可能有一种更参数化的方法来使用折刀残差计算置信区间。此代码来自S+ help page for smooth.spline
fit <- smooth.spline(data$x, data$y) # smooth.spline fit
res <- (fit$yin - fit$y)/(1-fit$lev) # jackknife residuals
sigma <- sqrt(var(res)) # estimate sd
upper <- fit$y + 2.0*sigma*sqrt(fit$lev) # upper 95% conf. band
lower <- fit$y - 2.0*sigma*sqrt(fit$lev) # lower 95% conf. band
matplot(fit$x, cbind(upper, fit$y, lower), type="plp", pch=".")
结果是
就 gam
置信区间而言,如果您阅读 print.gam
帮助文件,就会发现有一个 se=
参数默认为 TRUE
并且文档说
when TRUE (default) upper and lower lines are added to the 1-d plots at 2 standard errors above and below the estimate of the smooth being plotted while for 2-d plots, surfaces at +1 and -1 standard errors are contoured and overlayed on the contour plot for the estimate. If a positive number is supplied then this number is multiplied by the standard errors when calculating standard error curves or surfaces. See also shade, below.
所以可以通过调整这个参数来调整置信区间。 (这将在 print()
调用中进行。)
关于r - 如何获得 smooth.spline 的置信区间?,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/23852505/
我正在为我的雇主编写脚本,以从他们自己的站点获取某些数据。出于一长串原因,我需要从网站上获取数据,如图所示。我发现,其中一些数据是通过 js 调用检索的... 回想起来,我应该选择 Mechanize
我正在使用 python 和 cryptography.io 来签署和验证消息。我可以通过以下方式获得签名的 DER 编码字节表示: cryptography_priv_key.sign(messag
关闭。这个问题不满足Stack Overflow guidelines .它目前不接受答案。 想改善这个问题吗?更新问题,使其成为 on-topic对于堆栈溢出。 6年前关闭。 Improve thi
是否可以区分 ECDF?以下面得到的为例。 set.seed(1) a <- sort(rnorm(100)) b <- ecdf(a) plot(b) 我想对 b 求导以获得它的概率密度函数 (PD
我找到了如何从 navigator.mimeTypes 获取 mimetypes: function GetMimeTypes() { var message = ""; var mi
我在表单中使用单选按钮来隐藏/显示联系人表单中的成员 ID 字段。问题是,当 javascript 更改 html 中包含的隐藏 id 字段(该字段设置为“无”值)时,该字段将不再通过 post 可用
我正在做单元测试。我必须测试所有可能的if..else情况。但是在此if语句中: int32_t i32Res = snprintf(buffer, len, "The%d_String_%d", 0
我有一个 Facebook 应用程序,我想从中获取“喜欢”的总数。我想知道这是否可能。 其中 ID 是应用程序的 ID,ACCESS_TOKEN 是我尝试过的应用程序的当前访问 token : gra
如果我有多个计算实例尝试同时获取同一个 blob 的租约,则似乎经常会成功。我的印象是,一旦租约发出(并因此被客户获得),就不可能同时发出另一个租约? 我希望情况确实如此,我一直在 Azure 中使用
这是我的索引 POST /blogs/1 { "name" : "learn java", "popularity" : 100 } POST /blogs/2 { "name" : "l
我正在将 Symfony2 与 FOSUserBundle 一起使用。我需要为用户获得最高角色。 role_hierarchy: ROLE_CONTRIBUTOR: ROLE_USER
我正在向服务器发送基于 REST 的请求。我希望尽快得到答复,并希望了解可以进行的各种优化。 一种方法当然是在线程中并行发送这些请求。还有哪些其他选项可用于优化此功能? 在服务器上,可以添加哪些配置?
这可能是某种重复的问题,但我似乎找不到合适的解决方案。我正在使用 git4idea.history.GitHistoryUtils.history() 获取提交列表。如果 checkout 其中一个较
我正在做一个程序,可以输入每周的工资和那一周的总工作时间。它应该以小时工资率显示答案。但是我无法显示正确的“centavos/2 decimal places”公式并且它不想使用 float % fl
已结束。此问题正在寻求书籍、工具、软件库等的推荐。它不满足Stack Overflow guidelines 。目前不接受答案。 我们不允许提出寻求书籍、工具、软件库等推荐的问题。您可以编辑问题,以便
我已经尝试了 mContext.getMainLooper() 和 Looper.getMainLooper()。两者都返回相同的结果,但我想知道哪种方法正确? 我还从 Android 开发人员链接中
我有一个“affiliates”表,其中包含“user”和“referredBy”列。 给定一个用户,我希望获得该用户推荐的所有“n 级”玩家。对于 n=1,我们只关心您直接推荐的玩家数量: SELE
我在 PostgreSQL 9.5 数据库中有两个表: project - id - name task - id - project_id - name - updated_
请帮助我怎样才能得到我预期的结果,在此先感谢并抱歉我的英语不好。 PHP: $dog = implode(',', $data['dogbreed']); $query .= "AND `do
我有 let impulse = CGVectorMake(CGFloat(Constants.impulse), 0) 如何在不创建另一个 CGVector 的情况下得到它的负值? 我正在考虑在 C
我是一名优秀的程序员,十分优秀!