[R] confidence band
Hello,
I made a function, which calculates the confidence interval and the
prediction interval, but if I want to plot it, then it plots only the
regressionline.
Maybe somebody can help me:
conf.band - function(x,y) {
res.lsfit - lsfit(x,y)
xi - seq(from=40, to=160, length=200)
n - length(x)
MSE - sqrt(sum(res.lsfit$residuals^2)/(n-2))
# confidence interval
band - sqrt(2*qf(0.95,2,n-2)) * MSE *
sqrt(1/n+(xi-length(x))^2/(var(x)*(n-1)))
uiv - res.lsfit$coef[1] + res.lsfit$coef[2]*xi - band
oiv - res.lsfit$coef[1] + res.lsfit$coef[2]*xi + band
# prediction interval
band - sqrt(2*qf(0.95,2,n-2)) * MSE * sqrt(1+1/n
+(xi-mean(x))^2/(var(x)*(n-1)))
uip - res.lsfit$coef[1] + res.lsfit$coef[2]*xi - band
oip - res.lsfit$coef[1] + res.lsfit$coef[2]*xi + band
# creating the graphik
plot(x, y, xlab=colnames(x), ylab=colnames(y), pch=19)
abline(res.lsfit, col=1)
matlines(xi, cbind(uiv,oiv), col=3, lty=2, lwd=2)
matlines(xi, cbind(uiv,oip), col=2, lty=3, lwd=2)
}
--
Mit freundlichen Grüssen / Best Regards
Soare Marcian-Alin
[[alternative HTML version deleted]]
__
R-help@stat.math.ethz.ch 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.
RE: [R] Confidence Band for empirical distribution function
Dear Kjetil,
As I already mentioned, it appears that there isn't a function
available calculating the quantiles directly (at least, it doesn't appear
in the C source of ctest). So as I already suggested, uniroot (or a similar
C routine which calls the corresponding C code directly) is probably the
best you can do (apart from writing it completely yourself).
I didn't program this using uniroot, but I'd certainly try the following
for speed-up:
- For symmetry reasons, you only need to compute half of the quantiles.
- The quantiles depend smoothly on the probabilities (of your reference
distribution). Therefore, calculating only a few for probabilities
between 0 and 0.5, and using (e.g. linear) interpolation should be
satisfying.
I am sorry not be of more help.
HTH anyway
Thomas
-Original Message-
From: kjetil brinchmann halvorsen [mailto:[EMAIL PROTECTED]
Sent: 23 July 2003 15:46
To: Hotz, T.
Subject: RE: [R] Confidence Band for empirical distribution function
On 22 Jul 2003 at 11:37, Hotz, T. wrote:
Dear Leif,
If you look at the definition of ks.test, you'll find the lines
pkstwo - function(x, tol = 1e-06) {
if (is.numeric(x))
x - as.vector(x)
else stop(Argument x must be numeric)
p - rep(0, length(x))
p[is.na(x)] - NA
IND - which(!is.na(x) (x 0))
if (length(IND) 0) {
p[IND] - .C(pkstwo, as.integer(length(x)), p =
as.double(x[IND]),
as.double(tol), PACKAGE = ctest)$p
}
return(p)
}
which calls C code to calculate the p-values given the test
statistic.
You'll find explanations on what this function does in the
original C file
src/library/ctest/src/ks.c
I haven't tried that but I assume that you could use it to
calculate p-values
given the test-statistics yourself.
That could certainly be done, but what was asked for is the inverse,
which can be calculated, using for instance uniroot(). I tried that,
but it is to slow, .C() will be called repeatedly in a loop. For me
it took several minutes.
Kjetil Halvorsen
Please also note that ks.test() returns the p-value as well.
If you need quantiles, I assume you need to invert the cdf yourself,
e.g. using uniroot().
HTH
Thomas
---
Thomas Hotz
Research Associate in Medical Statistics
University of Leicester
United Kingdom
Department of Epidemiology and Public Health
22-28 Princess Road West
Leicester
LE1 6TP
Tel +44 116 252-5410
Fax +44 116 252-5423
Division of Medicine for the Elderly
Department of Medicine
The Glenfield Hospital
Leicester
LE3 9QP
Tel +44 116 256-3643
Fax +44 116 232-2976
-Original Message-
From: Leif.Boysen [mailto:[EMAIL PROTECTED]
Sent: 21 July 2003 14:42
To: [EMAIL PROTECTED]
Subject: [R] Confidence Band for empirical distribution function
Hi,
I was trying to draw an empirical distribution function
with uniform
confidence bands. So I tried to find a way to calculate
values of the
Kolmogorov-Smirnov Distribution but failed.
I guess it must be hidden somewhere (since the ks-test is
implemented),
but I was unable to find it.
Is there any way to do this?
Thanks
Leif Boysen
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
RE: [R] Confidence Band for empirical distribution function
On Thu, 24 Jul 2003, Hotz, T. wrote: Dear Kjetil, As I already mentioned, it appears that there isn't a function available calculating the quantiles directly (at least, it doesn't appear in the C source of ctest). So as I already suggested, uniroot (or a similar C routine which calls the corresponding C code directly) is probably the best you can do (apart from writing it completely yourself). I didn't program this using uniroot, but I'd certainly try the following for speed-up: - For symmetry reasons, you only need to compute half of the quantiles. - The quantiles depend smoothly on the probabilities (of your reference distribution). Therefore, calculating only a few for probabilities between 0 and 0.5, and using (e.g. linear) interpolation should be satisfying. There is an example of this in my package HyperbolicDist which has just appeared on CRAN. It is a little more sophisticated in that it fits a spline in preference to linear interpolation, before using uniroot. Look at qhyperb if this is of interest. David Scott _ David Scott Department of Statistics, Tamaki Campus The University of Auckland, PB 92019 AucklandNEW ZEALAND Phone: +64 9 373 7599 ext 86830 Fax: +64 9 373 7000 Email: [EMAIL PROTECTED] Graduate Officer, Department of Statistics Webmaster, New Zealand Statistical Association: http://www.stat.auckland.ac.nz/nzsa/ __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
RE: [R] Confidence Band for empirical distribution function
Dear Leif,
If you look at the definition of ks.test, you'll find the lines
pkstwo - function(x, tol = 1e-06) {
if (is.numeric(x))
x - as.vector(x)
else stop(Argument x must be numeric)
p - rep(0, length(x))
p[is.na(x)] - NA
IND - which(!is.na(x) (x 0))
if (length(IND) 0) {
p[IND] - .C(pkstwo, as.integer(length(x)), p = as.double(x[IND]),
as.double(tol), PACKAGE = ctest)$p
}
return(p)
}
which calls C code to calculate the p-values given the test statistic.
You'll find explanations on what this function does in the original C file
src/library/ctest/src/ks.c
I haven't tried that but I assume that you could use it to calculate p-values
given the test-statistics yourself.
Please also note that ks.test() returns the p-value as well.
If you need quantiles, I assume you need to invert the cdf yourself,
e.g. using uniroot().
HTH
Thomas
---
Thomas Hotz
Research Associate in Medical Statistics
University of Leicester
United Kingdom
Department of Epidemiology and Public Health
22-28 Princess Road West
Leicester
LE1 6TP
Tel +44 116 252-5410
Fax +44 116 252-5423
Division of Medicine for the Elderly
Department of Medicine
The Glenfield Hospital
Leicester
LE3 9QP
Tel +44 116 256-3643
Fax +44 116 232-2976
-Original Message-
From: Leif.Boysen [mailto:[EMAIL PROTECTED]
Sent: 21 July 2003 14:42
To: [EMAIL PROTECTED]
Subject: [R] Confidence Band for empirical distribution function
Hi,
I was trying to draw an empirical distribution function with uniform
confidence bands. So I tried to find a way to calculate values of the
Kolmogorov-Smirnov Distribution but failed.
I guess it must be hidden somewhere (since the ks-test is
implemented),
but I was unable to find it.
Is there any way to do this?
Thanks
Leif Boysen
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
[R] Confidence Band for empirical distribution function
Hi, I was trying to draw an empirical distribution function with uniform confidence bands. So I tried to find a way to calculate values of the Kolmogorov-Smirnov Distribution but failed. I guess it must be hidden somewhere (since the ks-test is implemented), but I was unable to find it. Is there any way to do this? Thanks Leif Boysen __ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
