[R] confidence intervals of proportions from complex surveys

2007-09-05 Thread Dirk Enzmann
This is partly an R and partly a general statistics question.

I'm trying to get confidence intervals of proportions (sometimes for 
subgroups) estimated from complex survey data. Because a function like 
prop.test() does not exist for the survey package I tried the following:

1) Define a survey object (PSU of clustered sample, population weights);
2) Use svyglm() of the package survey to estimate a binary logistic 
regression (family='binomial'): For the confidence interval of a single 
proportion regress the binary dependent variable on a constant (1), for 
confidence intervals of that variable for subgroups regress this 
variable on the groups (factor) variable;
3) Use predict() to obtain estimated logits and the respective standard 
errors (mod.dat specifiying either the constant or the subgroups):

pred=predict(model,mod.dat,type='link',se.fit=T)

and apply the following to obtain the proportion with its confidence 
intervals (for example, for conf.level=.95):

lo.e = pred[1:length(pred)]-qnorm((1+conf.level)/2)*SE(pred)
hi.e = pred[1:length(pred)]+qnorm((1+conf.level)/2)*SE(pred)
prop = 1/(1+exp(-pred[1:length(pred)]))
lo = 1/(1+exp(-lo.e))
hi = 1/(1+exp(-hi.e))

I think that in that way I get CI's based on asymptotic normality - 
either for a single proportion or split up into subgroups.

Question: Is this a correct or a defensible procedure? Or should I use a 
different approach? Note that this approach should also allow to 
estimate CI's for proportions of subgroups taking into account the 
complex survey design.

TIA,
Dirk


R version 2.5.1 Patched (2007-08-10 r42469)
i386-pc-mingw32

__
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.


[R] Confidence intervals for ccf()

2007-08-27 Thread Gustaf Rydevik
Hello,

This is not a purely R-question, but perhaps someone can help me anyway.

I am trying to estimate the correlation between two time series (which
are both basically different types of  measurements of the same
phenomena), using both cor.test() (with pearson as method) and ccf().

Now, cor.test gives a confidence interval for the pearson correlation,
while ccf does not. I've tried to use bootstrap methods to get
confidence interval for the ccf function, but no luck. It is a bit
tricky, since the time series are non-stationary, and so I'm not sure
how to go about to generate the bootstrap-sample.

Does anyone have any ideas on how to do this, i.e get a confidence
interval for the ccf at different time lags?

Many thanks in advance,

Gustaf

-- 
Gustaf Rydevik, M.Sci.
tel: +46(0)703 051 451
address:Essingetorget 40,112 66 Stockholm, SE
skype:gustaf_rydevik

__
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.


[R] confidence intervals for multinomial

2007-07-18 Thread Drescher, Michael (MNR)
Hi All,

I want to test an H0 hypothesis about the proportions of observed counts
in k classes. I know that I can do this with the chisq.test.

However, besides of the overall acceptance or rejection of the H0, I
would like to know which of the k classes cause(s) rejection and I would
like to know the observation-based confidence envelopes for the
proportions for the k classes.

My quick-and-dirty approach thus far is to do an initial chisq.test on
the original k classes and then to lump data into two classes (=one of
the original classes and all other original classes lumped into one new
class) and do a binom.test. I interpret the result of the binom.test as
indicating whether the current class might be the reason for the
rejection of the overall H0. Additionally, it gives me a confidence
envelope for this class.

This approach seems fairly straightforward, but I just do not feel
totally comfortable with it. I would feel so much better if there was
something like a multinom.test, but to my knowledge there is none.

Do you have any suggestions what I could rather do? For instance, I
might follow a Monte Carlo-like approach:
I simulate proportions for the k classes based on the proportions of
observed counts with rmultinom. After exclusion of the most extreme
values I construct my confidence envelope based on the remaining
simulated proportions. Based on whether the hypothesized proportions
fall into the observation-based confidence envelopes, I accept or
reject.

Do you think that either of these approaches is better or would you
suggest doing something totally different?

All comments and suggestions are highly appreciated.

Kind regards, Michael

PS: I guess my request parallels that of Matthias Schmidt from Apr 5,
2004, that was answered by Brian Ripley ...


Michael Drescher
Ontario Forest Research Institute
Ontario Ministry of Natural Resources
1235 Queen St East
Sault Ste Marie, ON, P6A 2E3
Tel: (705) 946-7406
Fax: (705) 946-2030

__
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 intervals on multiple comparisons

2007-05-15 Thread Cody_Hamilton

Enrico,

prop.test is for testing proportions two at a time.  If you want to test
for differences between 4 proportions simultaneously (rather than two at a
time), try a logistic regression model (from which you can get confidence
intervals for each of your groups).

Cody Hamilton, PhD
Staff Biostatistician
Edwards Lifesciences


   
 Salvatore Enrico 
 Indiogine
 [EMAIL PROTECTED]  To 
 .com R-help@stat.math.ethz.ch
 Sent by:   cc 
 [EMAIL PROTECTED] 
 at.math.ethz.ch   Subject 
   [R] confidence intervals on 
   multiple comparisons
 05/13/2007 10:51  
 AM
   
   
   
   




Greetings!

I am using prop.test to compare 4 proportions to find out whether they
are equal.  According to the help function you can not have confidence
intervals if you compare more than 2 proportions.

I need to find an effect size or confidence interval for these proportions.

Any suggestions?

Enrico

--
Enrico Indiogine

Mathematics Education
Texas AM University
[EMAIL PROTECTED]

__
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.

__
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 intervals on multiple comparisons

2007-05-15 Thread Bert Gunter
 

 -Original Message-
 From: [EMAIL PROTECTED] 
 [mailto:[EMAIL PROTECTED] On Behalf Of 
 [EMAIL PROTECTED]
 Sent: Tuesday, May 15, 2007 12:52 PM
 To: Salvatore Enrico Indiogine
 Cc: R-help@stat.math.ethz.ch; [EMAIL PROTECTED]
 Subject: Re: [R] confidence intervals on multiple comparisons
 
 
 Enrico,
 
 prop.test is for testing proportions two at a time.  If you 
 want to test
 for differences between 4 proportions simultaneously (rather 
 than two at a
 time), try a logistic regression model (from which you can 
 get confidence
 intervals for each of your groups).
 
 Cody Hamilton, PhD
 Staff Biostatistician
 Edwards Lifesciences
 

Yes, but beware: in the default contr.treatment coding for contrasts, you
get estimates and confidence intervals for the first group and for the
**differences** between the first group and the others. As you said, it's
easy to get what you want from this, but you must pay attention to the
details here. 

Bert Gunter
Genentech Nonclinical statistics


   
  
  Salvatore Enrico
  
  Indiogine   
  
  [EMAIL PROTECTED]
   To 
  .com 
 R-help@stat.math.ethz.ch
  Sent by: 
   cc 
  [EMAIL PROTECTED]
  
  at.math.ethz.ch  
  Subject 
[R] confidence 
 intervals on 
multiple comparisons   
  
  05/13/2007 10:51 
  
  AM   
  
   
  
   
  
   
  
   
  
 
 
 
 
 Greetings!
 
 I am using prop.test to compare 4 proportions to find out whether they
 are equal.  According to the help function you can not have confidence
 intervals if you compare more than 2 proportions.
 
 I need to find an effect size or confidence interval for 
 these proportions.
 
 Any suggestions?
 
 Enrico
 
 --
 Enrico Indiogine
 
 Mathematics Education
 Texas AM University
 [EMAIL PROTECTED]
 
 __
 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.
 
 __
 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.


__
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.


[R] Confidence-Intervals.... help...

2007-04-08 Thread Jochen.F

Hi...

I have to use R to find out the 90% confidence-interval for the sensitivity
and specificity of the following diagnostic test:

A particular diagnostic test for multiple sclerosis was conducted on 20 MS
patients and 20 healthy subjects, 6 MS patients were classified as healthy
and 8 healthy subjects were classified as suffering from the MS.

Furthermore, I need to find the number of MS patients required for a
sensitivity of 1%...

Is there a simple R-command which can do that for me?

I am completely new to R...

Help please!

Jochen
-- 
View this message in context: 
http://www.nabble.com/Confidence-Intervals-help...-tf3544217.html#a9894014
Sent from the R help mailing list archive at Nabble.com.

__
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-Intervals.... help...

2007-04-08 Thread Sarah Goslee
Hm... sounds like a homework problem to me...

Maybe start by figuring out how to do it without R -
what's the approach, and how would you calculate it?
Then search R help for the possible key words you
came up with.

Sarah

On 4/8/07, Jochen.F [EMAIL PROTECTED] wrote:

 Hi...

 I have to use R to find out the 90% confidence-interval for the sensitivity
 and specificity of the following diagnostic test:

 A particular diagnostic test for multiple sclerosis was conducted on 20 MS
 patients and 20 healthy subjects, 6 MS patients were classified as healthy
 and 8 healthy subjects were classified as suffering from the MS.

 Furthermore, I need to find the number of MS patients required for a
 sensitivity of 1%...

 Is there a simple R-command which can do that for me?

 I am completely new to R...

 Help please!

 Jochen
 --
-- 
Sarah Goslee
http://www.functionaldiversity.org

__
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.


[R] confidence intervals

2007-02-22 Thread Bart Joosen
Hi,

I'm having trouble with the confidence interval of the nls function.
I did my home work, and searched acros the support list until I came up with 
following solution of Peter Dalgaard:

example(predict.nls) 
se.fit - sqrt(apply(attr(predict(fm,list(Time = 
tt)),gradient),1, 
function(x) sum(vcov(fm)*outer(x,x
matplot(tt, predict(fm,list(Time = tt))+
   outer(se.fit,qnorm(c(.5, .025,.975))),type=l)
points(demand ~ Time, data = BOD) 
One slight issue is that it doesn't work if newdata is omitted, 
but then you can easily get the gradient from fm$m$gradient() 

I tried this with my own data:
Data - data.frame(Temp=rep(c(25,40),each=3), Mnd = c(1:3),Degr = 
c(0.057,0.077,0.108,0.148,0.198,0.223))
model - nls(Degr~exp((A/(Temp)+log(Mnd))*B),Data,start=list(A=-10,B=1))
Months - c(1,2,3,6,9,12,24,48)
se.fit - sqrt(apply(attr(predict(model,list(Temp = 
25,Mnd=Months)),gradient),1, function(x) sum(vcov(fm)*outer(x,x

But unfortunately I get an error ( Error in apply(attr(predict(model, list(Temp 
= 25, Mnd = Months)), gradient),  : 
dim(X) must have a positive length)

When I try using the gradient of the model instead of using the new data then 
there is no problem:
se.fit - sqrt(apply(model$m$gradient(),1, function(x) 
sum(vcov(model)*outer(x,x
matplot(Data$Mnd, predict(model,list(Temp = 
Data$Temp,Mnd=Data$Mnd))+outer(se.fit,qnorm(c(.5, 025,.975))),type=l)

But how about calculating confidence intervals of new data? How do I get an 
gradient for these values?

I'm using Windows XP, R 2.4.1.


Thanks

Bart
[[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 intervals of quantiles

2007-02-09 Thread mikewhite . diu
Many thanks for all the contributions to this problem.
As inferred by Ted Harding, I was after a distribution-free CIs as a lot
of
the data I use is not normally distributed.

The method provided by Ted for calculating exact CIs gave good results with
the limits almost symmetric about the quantile.  Even for my smaller data
set with only 115 samples and a very skewed distribution the results clearly
showed the increase in range and asymmetry of the confidence limits.

The bootstrap method provided by Dimitris also works reasonably well but
is
slower and the ranges for the CIs are sometimes very asymmetric and in one
case did not actually encompass the quantile. The method would not work at
all with the skewed distribution of 115 samples until I reduced the quantile
range from 0.1 and 0.9 to 0.2 and 0.8  However, as Dimitris warned, you have
to be careful with this method for extreme quantiles and small samples.

I was also sent a method by Søren Merser, but this was incomplete.  This
method was written by Scott Chasalow and the full code can be found at
http://www.dpw.wau.nl/pv/pub/chasalow/S/win/ci.quantile/
The code was written for S-Plus but it worked ok for me in R.  This method
actually gives several ranges about the quantile, each with about the same
level of confidence and the level of confidence is also in the output. As
with Ted Harding's method, these may not exactly match the desired
confidence level.  There is an option to select the shortest range but it
would be easy enough to add code to give the most symmetric range.

As a chemist I am not able to comment on the statistical pros and cons of
the methods but they are certainly very helpful for my purposes.

Many thanks
Mike White

- Original Message - 
From: Dimitris Rizopoulos [EMAIL PROTECTED]
To: Mike White [EMAIL PROTECTED]
Cc: R-help@stat.math.ethz.ch
Sent: Monday, February 05, 2007 2:43 PM
Subject: Re: [R] Confidence intervals of quantiles


 you could use the Bootstrap method, using package 'boot', e.g.,

 library(boot)

 f.quantile - function(x, ind, ...){
 quantile(x[ind], ...)
 }

 ###

 x - rgamma(750, 2)
 quant.boot - boot(x, f.quantile, R = 1000, probs = c(0.025, 0.25,
 0.5, 0.75, 0.975))
 lapply(1:5, function(i) boot.ci(quant.boot, c(0.90, 0.95), type =
 c(perc, bca), index = i))

 y - rgamma(150, 2)
 quant.boot - boot(y, f.quantile, R = 1000, probs = c(0.025, 0.25,
 0.5, 0.75, 0.975))
 lapply(1:5, function(i) boot.ci(quant.boot, c(0.90, 0.95), type =
 c(perc, bca), index = i))


 However, you should be a little bit careful with Bootstrap if you wish
 to obtain CIs for extreme quantiles in small samples.

 I hope it helps.

 Best,
 Dimitris

 
 Dimitris Rizopoulos
 Ph.D. Student
 Biostatistical Centre
 School of Public Health
 Catholic University of Leuven

 Address: Kapucijnenvoer 35, Leuven, Belgium
 Tel: +32/(0)16/336899
 Fax: +32/(0)16/337015
 Web: http://med.kuleuven.be/biostat/
  http://www.student.kuleuven.be/~m0390867/dimitris.htm


 - Original Message - 
 From: Mike White [EMAIL PROTECTED]
 To: R-help@stat.math.ethz.ch
 Sent: Monday, February 05, 2007 2:47 PM
 Subject: [R] Confidence intervals of quantiles


  Can anyone please tell me if there is a function to calculate
  confidence
  intervals for the results of the quantile function.
  Some of my data is normally distributed but some is also a squewed
  distribution or a capped normal distribution. Some of the data sets
  contain
  about 700 values whereas others are smaller with about 100-150
  values, so I
  would like to see how the confidence intervals change for the
  different
  distributions and different data sizes.
 
  Thanks
  Mike White
 
  __
  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.
 


 Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm




___

Tiscali Broadband only 9.99 a month for your first 3 months!
http://www.tiscali.co.uk/products/broadband/

__
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.


[R] Confidence intervals of quantiles

2007-02-05 Thread Mike White
Can anyone please tell me if there is a function to calculate confidence
intervals for the results of the quantile function.
Some of my data is normally distributed but some is also a squewed
distribution or a capped normal distribution. Some of the data sets contain
about 700 values whereas others are smaller with about 100-150 values, so I
would like to see how the confidence intervals change for the different
distributions and different data sizes.

Thanks
Mike White

__
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 intervals of quantiles

2007-02-05 Thread Dimitris Rizopoulos
you could use the Bootstrap method, using package 'boot', e.g.,

library(boot)

f.quantile - function(x, ind, ...){
quantile(x[ind], ...)
}

###

x - rgamma(750, 2)
quant.boot - boot(x, f.quantile, R = 1000, probs = c(0.025, 0.25, 
0.5, 0.75, 0.975))
lapply(1:5, function(i) boot.ci(quant.boot, c(0.90, 0.95), type = 
c(perc, bca), index = i))

y - rgamma(150, 2)
quant.boot - boot(y, f.quantile, R = 1000, probs = c(0.025, 0.25, 
0.5, 0.75, 0.975))
lapply(1:5, function(i) boot.ci(quant.boot, c(0.90, 0.95), type = 
c(perc, bca), index = i))


However, you should be a little bit careful with Bootstrap if you wish 
to obtain CIs for extreme quantiles in small samples.

I hope it helps.

Best,
Dimitris


Dimitris Rizopoulos
Ph.D. Student
Biostatistical Centre
School of Public Health
Catholic University of Leuven

Address: Kapucijnenvoer 35, Leuven, Belgium
Tel: +32/(0)16/336899
Fax: +32/(0)16/337015
Web: http://med.kuleuven.be/biostat/
 http://www.student.kuleuven.be/~m0390867/dimitris.htm


- Original Message - 
From: Mike White [EMAIL PROTECTED]
To: R-help@stat.math.ethz.ch
Sent: Monday, February 05, 2007 2:47 PM
Subject: [R] Confidence intervals of quantiles


 Can anyone please tell me if there is a function to calculate 
 confidence
 intervals for the results of the quantile function.
 Some of my data is normally distributed but some is also a squewed
 distribution or a capped normal distribution. Some of the data sets 
 contain
 about 700 values whereas others are smaller with about 100-150 
 values, so I
 would like to see how the confidence intervals change for the 
 different
 distributions and different data sizes.

 Thanks
 Mike White

 __
 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.
 


Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm

__
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 intervals of quantiles

2007-02-05 Thread Ted Harding
On 05-Feb-07 Mike White wrote:
 Can anyone please tell me if there is a function to calculate
 confidence intervals for the results of the quantile function.
 Some of my data is normally distributed but some is also a
 squewed distribution or a capped normal distribution. Some of
 the data sets contain about 700 values whereas others are smaller
 with about 100-150 values, so I would like to see how the confidence
 intervals change for the different distributions and different data
 sizes.

As well as the bootstrap suggestion (which may do what you want)
I'd like to suggest calcualting exact CIs for the distribution
quantiles. I don't know of an R function which does this, though
the principle is straightforward enough (I once implemented it
for octave/matlab).

I'm assuming you're after a truly distribution-free CI (as your
query suggests). In that case, the sample quantiles provide the
CI limits for the distribution quantiles, with the proviso that
you will not in general get an exact match for the confidence
level P that you desire for your confidence interval, so you
need to use P as a lower bound for the confidence level.

I'll illustrate with an equi-tailed  95% CI.

Notation:

x[r] (r = 1:n) is the r-th order statistic of a sample of n
values of a random variable X, drawn from a continuous distribution.
X[r] is the corresponding random variable.

X(p) is the p-th quantile of the distribution in question,
i.e. the value of X such that P(X = X[p]) = p. 0  p  1.

Objective: You want a 95% CI (X(p)L,X(p).R) for the quantile X(p)
for given p.

Principle:

  P( X[r] = X(p) ) = pbinom(r, n, p)  [*]

(i.e. the probability that you get at least r of the x's
below X(p)).

Now, for the upper limit X(p).R of the CI, you want the
smallest value of X(p) such that the probability [*] is not
less than 0.925 (i.e. 1 - (1-P)/2 = 1 - (1-0.95)/2 ).

This is achieved by X(p).R = x[s] where

  s = min(which(pbinom((0:n),n,p) = 0.025))

and, similalry, for the lower limit, X(p).L = x[r] where

  r = max(which(pbinom((0:n),n,p) = 0.025))

(Note that the which counts the binomial case r=0 as number 1)

If I've got my head around the above right, the following
function does the above job, and caould be fairly easily
modified for asymmetric confidence intervals (e.g. a 95%
confidence interval with 96% for inclusion below the upper
limit and 99% for inclusion above the lower limit).

q.CI-function(x,p,P){
# x is the sample, p (0p1) the quantile, P the confidence level
  x-sort(x)
  n-length(x)
  s - min(which(pbinom((0:n),n,p) = 1-(1-P)/2))
  r - max(which(pbinom((0:n),n,p) = (1-P)/2))
  c(x[r],x[s])
# x[r] is the lower limit, x[s] the upper limit, of the CI
}

Note that it gives a confidence level P at least equal to P,
not in general exactly equal to P.

I've tested it as follows (1 simulations with samples
of size 101 from a Normal distribution -- why not, after
all? No loss of generality!):

p-0.5; Xq-qnorm(p); P-0.95
N-1; incls-numeric(N)
for(i in (1:N)){
  x-rnorm(101)
  CI-q.CI(x,p,P)
  if((CI[1]=Xq)(CI[2]=Xq)){ incls[i]-1 }
}; sum(incls)/N
# [1] 0.9564
# [1] 0.9548

p-0.75; Xq-qnorm(p); P-0.95
[etc. ... ]
# [1] 0.9631
# [1] 0.9596

p-0.90; Xq-qnorm(p); P-0.95
[etc. ... ]
# [1] 0.9540
# [1] 0.9533

p-0.95; Xq-qnorm(p); P-0.95
[etc. ... ]
# [1] 0.9840
# [1] 0.9826

p-0.99; Xq-qnorm(p); P-0.95
[etc. ... ]
Error in if ((CI[1] = Xq)  (CI[2] = Xq)) { : 
missing value where TRUE/FALSE needed

showing that for extreme values of p relative to the
sample size, trouble occurs (as is to be expected)!

The solution is to test for NA in CI[1] and/or CI[2],
so I modified my test routine as follows, to give notional
lower confidence limit -Inf if CI[1] is NA, and/or upper
confidence limit +Inf if CI[2] is NA (of course this
could be done in the function q.CI() itself, but perhaps
it offers more flexibility for the user to do something
else with it, if it is left as NA).

p-0.99; Xq-qnorm(p); P-0.95
N-1; incls-numeric(N)
for(i in (1:N)){
  x-rnorm(101)
  CI-q.CI(x,p,P)
  if(is.na(CI[1])){ CI[1]-(-Inf) }
  if(is.na(CI[2])){ CI[2]-(+Inf) }
  if((CI[1]=Xq)(CI[2]=Xq)){ incls[i]-1 }
}; sum(incls)/N
# [1] 0.9841
# [1] 0.9821

Comments welcome!!!

Ted.




E-Mail: (Ted Harding) [EMAIL PROTECTED]
Fax-to-email: +44 (0)870 094 0861
Date: 05-Feb-07   Time: 23:40:23
-- XFMail --

__
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 intervals for relative risk

2006-11-13 Thread Terry Therneau

Wolfgang,
  It is common to handle relative risk problems using Poisson regression.
  In your example you have 8 events out of 508 tries, and 0/500 in the second
data set.

 tdata - data.frame(y=c(8,0), n=c(508,500), group=1:0)
 fit - glm(y ~ group + offset(log(n)), data=tdata, family=poisson)

  Because of the zero, the standard beta/se(beta) confidence intervals don't
work.  In fact, the actual MLE for the ratio is infinity, which the above
glm decides is exp(33.85) = 5* 10^14 --- which is close enough to infinity
for me.  What you want is the lower limit of the interval, which you can
find with a profile likelihood.  That is, look at the curve of deviance vs
beta, draw a horizontal line 3.84 units down from the top (chisq on 1 df),
and where it itersects the curve is your confidence limit.

  For the above, since it is 2 groups with 2 coefficients, the deviance of
the full model happens to be 0.  Your approximation gave a lower limit of
.97 which is about exp(0), so I'll guess a solution between -2 and 2.

 xx - seq(-2, 2, length=21)
 for (i in 1:21) {
fit - glm(y ~ offset(group*xx[i] + log(n)), poisson, tdata)
print(c(xx[i], fit$deviance))
}

[1]  -2.0  33.80736
[1]  -1.8  31.02994
[1]  -1.6  28.33138
[1]  -1.4  25.72327
[1]  -1.2  23.21769
[1]  -1.0  20.82691
[1]  -0.8  18.56281
[1]  -0.6  16.43629
[1]  -0.4  14.45668
[1]  -0.2  12.63108
[1]  0.0 10.96387
[1] 0.2 9.45639
[1] 0.40 8.106805
[1] 0.60 6.910274
[1] 0.80 5.859303
[1] 1.00 4.944278
[1] 1.20 4.154088
[1] 1.40 3.476757
[1] 1.6 2.90003
[1] 1.8 2.41186
[1] 2.00 2.000785

 So the exact lower limit is somewhere between exp(1.4) and exp(1.2), which
is around 3.6.  One can easily refine this by tucking the fit into an
iterative search, or just resetting xx to a new range.

  The offset(log(n)) part of the model, BTW, is a well known trick in poisson
models.  It has to do with the fact that the likelihood is in terms of the
number of events, but we want coefficients in terms of rates, and E(y) = rate*n.
Adding an offset of xx[i] * group fits a model with the group effect fixed at
xx, but allowing the intercept to vary.  
  For more conceptual depth, you could look up 'rate regression' or 
'standardized mortality ratio' in the Encyclopedia of Biostatistics --- it is
in the latter computation that this approach is quite common.  The idea of
adding a fraction is found in Anscombe(1949), Transformations of Poisson, 
binomial and negative-binomial data.  Biometrika, vol 35, p246-254, but for
each single estimate not for the ratio.

Terry Therneau
Mayo Clinic

__
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.


[R] Confidence intervals for Brier score

2006-11-01 Thread Inman, Brant A. M.D.

I am using the ipred library to calculate the censored Brier score for a
Cox proportional hazards model.  I would like to know if anyone has
developed a method of calculating confidence intervals for the various
forms of the Brier score that are used in the analysis of
survival/censored data.  If so, are these implemented in S?

Thank you

Brant Inman

[[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.


[R] Confidence intervals on Lowess (eg SiZer from J.S.Marron)

2006-09-08 Thread Dalphin, Mark
Hi,

I have some very noisy, relatively sparse data; a biological response of
roughly
~8 subjects at ~8 times points). I've been following the data trend using a
lowess
line, over-plotted with several values of bandwidth, 'f - seq(0.3, 0.9,
by=0.1)'.
At this point, we have no models for these data.

I wonder if there is any way under R to assign some sort of confidence
interval to the
lowess line. For example, I have seen a method from the lab of J.S.Marron,
which
his group has implemented in the Matlab program, SiZer.
http://www.stat.unc.edu/faculty/marron/DataAnalyses/SiZer_Intro.html
and, the more interesting:

http://www.stat.unc.edu/faculty/marron/DataAnalyses/SiZer/SiZer_Basics.html

An implementation under R of SiZer would obviously answer this question;
does it exist?
Suggestions of alternative approaches would also be welcome.

My searchs of the R-maillist archive for 'sizer', 'Marron' and 'scale-space'
didn't return anything
that I recognized as useful here, but I am new to the idea of confidence in
a lowess line
so I may not be using the appropriate vocabulary.

Thank-you,
Mark

--
Mark Dalphin
Dept Comp Biol, M/S AW2/D3262
Amgen, Inc.
1201 Amgen Court W
Seattle, WA 98119
Phone: +1-206-265-7951

__
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 Intervals for Mixed Effects

2005-11-05 Thread Spencer Graves
  Are you familiar with intervals and lme in the nlme package? 
These are documented in Pinheiro and Bates (2000) Mixed-Effects Models 
in S and S-Plus (Springer).  I'm not familiar with the algorithm, but if 
it's different from Satterthwaite's method, I suspect that Prof. Bates 
had a good reason for choosing something different.  Since R is open 
source, you could read the code and modify it to use Sattherthwaite and 
compare the two side by side.  The nlme package includes a 
simulate.lme function, which you could use to compare the different 
methods.

  hope this helps.
  spencer graves

Michel Friesenhahn wrote:
 I'm fairly new to R and am wondering if anybody knows of R code to 
 calculate confidence intervals for parameters (fixed effects and variance 
 components) from mixed effects models based on Sattherthwaite's method? 
 I'm also interested in Satterthwaite-based confidence intervals for linear 
 combinations (mostly sums) of various variance components.
 
   [[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

-- 
Spencer Graves, PhD
Senior Development Engineer
PDF Solutions, Inc.
333 West San Carlos Street Suite 700
San Jose, CA 95110, USA

[EMAIL PROTECTED]
www.pdf.com http://www.pdf.com
Tel:  408-938-4420
Fax: 408-280-7915

__
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


[R] Confidence Intervals for Mixed Effects

2005-10-25 Thread Michel Friesenhahn
I'm fairly new to R and am wondering if anybody knows of R code to 
calculate confidence intervals for parameters (fixed effects and variance 
components) from mixed effects models based on Sattherthwaite's method? 
I'm also interested in Satterthwaite-based confidence intervals for linear 
combinations (mostly sums) of various variance components.

[[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


Re: [R] Confidence Intervals for Arbitrary Functions

2005-07-17 Thread Spencer Graves
  Before you spend a lot of time worrying about the distributions of 
ratios of normal variates, bootstrapping, etc., I suggest you make 
normal probability plots of your data AND OF THEIR LOGARITHMS, e.g., using

qqnorm(y, datax=TRUE)
qqnorm(y, datax=TRUE, log=x)

  With many electrical measurements from integrated circuits, these two 
plots will both likely be equally close to normality.  That's not true 
for leakage currents, for which I often see highly skewed images, unless 
I take logarithms, as:

qqnorm(exp(rnorm(99)), datax=TRUE)
qqnorm(exp(rnorm(99)), datax=TRUE, log=x)

  I routinely see distributions of this nature in estimates of Poisson 
defect rates in wafer fab as well as in leakage currents, distributions 
of income, etc.

  This is important, because the distribution of a ratio of normally 
distributed random variables has known pathologies, while the 
distribtuion of a ratio of lognormal variates is lognormal.  (With a 
ratio of normals, if the denominator has mean zero, the ratio follows 
the Cauchy disribution, also known as Student's t with one degree of 
freedom.  This distribution has infinite variance, and the mean is not 
even defined, being Inf-Inf.)

  In sum, you will likely get better answers in less time if you can 
find it politically acceptable in your professional efforts to work in 
decibels or logaritms, where ratios become differences.

  spencer graves

Gabor Grothendieck wrote:

 On 7/16/05, Jeff Newmiller [EMAIL PROTECTED] wrote:
 
I have a rather basic background in statistics, and am looking for
assistance in solving what I expect is a common type of problem.

I have measurements of physical processes, and mathematical models of
those processes that I want to feed the measurements into. A simple case
is using measurements of electric power entering and leaving a
power conversion device, sampled at regular intervals, and summed to
estimate energy in and out, and dividing the energy out by the energy in
to get an estimate of efficiency.  I know that power efficiency varies
with power level, but for this calculation I am interested in the
quantifying the overall efficiency rather than the instantaneous
efficiency.

If the energy quantities are treated as a normally-distributed random
variable (per measurement uncertainty), is there a package that simplifies
the determination of the probability distribution function for the
quotient of these values? Or, in the general sense, if I have a function
that computes a measure of interest, are such tools general enough to
handle this? (The goal being to determine a confidence interval for the
computed quantity.)

As an attempt to understand the issues, I have used SQL to generate
discrete sampled normal distributions, and then computed new abscissa
values using a function such as division and computing the joint
probability as the ordinate, and then re-partitioned the result into new
bins using GROUP BY.  This is general enough to handle non-normal
distributions as well, though I don't know how to quantify the numerical
stability/accuracy of this computational procedure. However, this is
pretty tedious... it seems like R ought to have some straightforward
solution to this problem, but I don't seem to know what search terms to
use.

 
 
 There is some discussion about the ratio of normals at:
 
http://www.pitt.edu/~wpilib/statfaq.html
 
 but you may just want to use bootstrapping:
 
   library(boot)
   library(simpleboot)
 
 __
 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

-- 
Spencer Graves, PhD
Senior Development Engineer
PDF Solutions, Inc.
333 West San Carlos Street Suite 700
San Jose, CA 95110, USA

[EMAIL PROTECTED]
www.pdf.com http://www.pdf.com
Tel:  408-938-4420
Fax: 408-280-7915

__
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


[R] Confidence Intervals for Arbitrary Functions

2005-07-16 Thread Jeff Newmiller
I have a rather basic background in statistics, and am looking for
assistance in solving what I expect is a common type of problem.

I have measurements of physical processes, and mathematical models of
those processes that I want to feed the measurements into. A simple case
is using measurements of electric power entering and leaving a
power conversion device, sampled at regular intervals, and summed to
estimate energy in and out, and dividing the energy out by the energy in
to get an estimate of efficiency.  I know that power efficiency varies
with power level, but for this calculation I am interested in the
quantifying the overall efficiency rather than the instantaneous
efficiency.

If the energy quantities are treated as a normally-distributed random
variable (per measurement uncertainty), is there a package that simplifies
the determination of the probability distribution function for the
quotient of these values? Or, in the general sense, if I have a function
that computes a measure of interest, are such tools general enough to
handle this? (The goal being to determine a confidence interval for the
computed quantity.)

As an attempt to understand the issues, I have used SQL to generate
discrete sampled normal distributions, and then computed new abscissa
values using a function such as division and computing the joint
probability as the ordinate, and then re-partitioned the result into new
bins using GROUP BY.  This is general enough to handle non-normal
distributions as well, though I don't know how to quantify the numerical
stability/accuracy of this computational procedure. However, this is
pretty tedious... it seems like R ought to have some straightforward
solution to this problem, but I don't seem to know what search terms to
use.

---
Jeff NewmillerThe .   .  Go Live...
DCN:[EMAIL PROTECTED]Basics: ##.#.   ##.#.  Live Go...
  Live:   OO#.. Dead: OO#..  Playing
Research Engineer (Solar/BatteriesO.O#.   #.O#.  with
/Software/Embedded Controllers)   .OO#.   .OO#.  rocks...1k

__
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


Re: [R] Confidence Intervals for Arbitrary Functions

2005-07-16 Thread Gabor Grothendieck
On 7/16/05, Jeff Newmiller [EMAIL PROTECTED] wrote:
 I have a rather basic background in statistics, and am looking for
 assistance in solving what I expect is a common type of problem.
 
 I have measurements of physical processes, and mathematical models of
 those processes that I want to feed the measurements into. A simple case
 is using measurements of electric power entering and leaving a
 power conversion device, sampled at regular intervals, and summed to
 estimate energy in and out, and dividing the energy out by the energy in
 to get an estimate of efficiency.  I know that power efficiency varies
 with power level, but for this calculation I am interested in the
 quantifying the overall efficiency rather than the instantaneous
 efficiency.
 
 If the energy quantities are treated as a normally-distributed random
 variable (per measurement uncertainty), is there a package that simplifies
 the determination of the probability distribution function for the
 quotient of these values? Or, in the general sense, if I have a function
 that computes a measure of interest, are such tools general enough to
 handle this? (The goal being to determine a confidence interval for the
 computed quantity.)
 
 As an attempt to understand the issues, I have used SQL to generate
 discrete sampled normal distributions, and then computed new abscissa
 values using a function such as division and computing the joint
 probability as the ordinate, and then re-partitioned the result into new
 bins using GROUP BY.  This is general enough to handle non-normal
 distributions as well, though I don't know how to quantify the numerical
 stability/accuracy of this computational procedure. However, this is
 pretty tedious... it seems like R ought to have some straightforward
 solution to this problem, but I don't seem to know what search terms to
 use.
 

There is some discussion about the ratio of normals at:

   http://www.pitt.edu/~wpilib/statfaq.html

but you may just want to use bootstrapping:

  library(boot)
  library(simpleboot)

__
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


[R] Confidence intervals for rates (dependent events)

2005-02-08 Thread Dirk Enzmann
I need advice or opinions for the following problem:
In a sample a part of the respondents has experienced victimizing 
events. Only a part of these events have been reported to the police. 
The rate of events reported to the police is

number of experienced events / number of reported events.
Because the events cannot be assumed to be independent (some victims 
have a high rate of vicitimization because they are more prone to become 
a victim) the construction of a confidence interval for the rate of 
events reported becomes difficult (to me).

Question 1: If nevertheless a use a binomial test to construct a 
confidence interval (say: 0 of 13 events reported,

binom.test(0,13,0/13)
CI = 0 to 24.7 %), is it correct that the width of this interval is a 
lower bound and thus a conservative estimate?

Question 2: (a) My intuition tells me that multilevel modeling could be 
a solution to obtain a correct confidence interval by treating the 
events and the events reported as the first level and the victims as the 
second. Is this correct and how should I specify this model?

(b) Alternatively, could I treat the events (and events reported?) as 
coming from a negative binomial distribution and use this for 
constructing a confidence interval? How can this be done technically, 
for example by using nb.glm?

Thanks in advance,
Dirk
*
Dr. Dirk Enzmann
Institute of Criminal Sciences
Dept. of Criminology
Schlueterstr. 28
D-20146 Hamburg
Germany
phone: +49-040-42838.7498 (office)
   +49-040-42838.4591 (Billon)
fax:   +49-040-42838.2344
email: [EMAIL PROTECTED]
www: 
http://www2.jura.uni-hamburg.de/instkrim/kriminologie/Mitarbeiter/Enzmann/Enzmann.html

__
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


RE: [R] Confidence Intervals from Bootstrap Replications

2004-12-19 Thread Ashraf Chaudhary
Thank you for the response. The suggestions are helpful. I would look into
the capabilities of the boot function. The solution in 8.3 is the last
resort that I was trying to avoid. Regards,
Ashraf

-Original Message-
From: Berton Gunter [mailto:[EMAIL PROTECTED] 
Sent: Friday, December 17, 2004 5:54 PM
To: 'Mohammad A. Chaudhary'; [EMAIL PROTECTED]
Subject: RE: [R] Confidence Intervals from Bootstrap Replications 


See , e.g. section 8.3 The two-sample problem of Efron and Tibshirani's AN
INTRODUCTION TO THE BOOTSTRAP. It makes it clear there why one just
bootstrap samples independently from the two separate samples.

The strata argument of boot() in the boot package allows one to do such
independent sampling and use the capabilities of that function.

-- Bert Gunter
Genentech Non-Clinical Statistics
South San Francisco, CA
 
The business of the statistician is to catalyze the scientific learning
process.  - George E. P. Box
 
 

 -Original Message-
 From: [EMAIL PROTECTED] 
 [mailto:[EMAIL PROTECTED] On Behalf Of 
 Mohammad A. Chaudhary
 Sent: Friday, December 17, 2004 2:39 PM
 To: [EMAIL PROTECTED]
 Subject: [R] Confidence Intervals from Bootstrap Replications 
 
 Hi All:
 
 I have to compute bootstrap confidence intervals, the statistic
 (incremental cost effectiveness ratio) is computed from two samples
 (intervention and control) of different sizes. All the bootstrap
 functions that I have seen use one dataset as argument. I may go ahead
 and get the desired number of bootstrap replications 
 separately. I would
 appreciate if you could point me to a source of a bootstrap 
 function (if
 available) that takes the B bootstrap replications and other 
 descriptive
 statistics and can get me the confidence intervals. Please 
 write me if I
 have not been clear in explaining my problem. Regards,
 
 Ashraf  
 
  
 
 ___
 
 M. Ashraf Chaudhary, Ph.D.
 
 Associate Scientist/Biostatistician
 
 Department of International Health
 
 Disease Prevention and Control Program
 
 Johns Hopkins University Bloomberg School of Public Health
 
 615 North Wolfe Street, Room W5506
 
 Baltimore MD 21205
 
  
 
 [EMAIL PROTECTED] 
 
 Phone: (410) 502-0741/Fax: (410) 502-6733
 
 http://faculty.jhsph.edu/?F=MohammadL=Chaudhary
 
  
 
  
 
 
   [[alternative HTML version deleted]]
 
 __
 [EMAIL PROTECTED] mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! 
 http://www.R-project.org/posting-guide.html


__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] Confidence Intervals from Bootstrap Replications

2004-12-17 Thread Mohammad A. Chaudhary
Hi All:

I have to compute bootstrap confidence intervals, the statistic
(incremental cost effectiveness ratio) is computed from two samples
(intervention and control) of different sizes. All the bootstrap
functions that I have seen use one dataset as argument. I may go ahead
and get the desired number of bootstrap replications separately. I would
appreciate if you could point me to a source of a bootstrap function (if
available) that takes the B bootstrap replications and other descriptive
statistics and can get me the confidence intervals. Please write me if I
have not been clear in explaining my problem. Regards,

Ashraf  

 

___

M. Ashraf Chaudhary, Ph.D.

Associate Scientist/Biostatistician

Department of International Health

Disease Prevention and Control Program

Johns Hopkins University Bloomberg School of Public Health

615 North Wolfe Street, Room W5506

Baltimore MD 21205

 

[EMAIL PROTECTED] 

Phone: (410) 502-0741/Fax: (410) 502-6733

http://faculty.jhsph.edu/?F=MohammadL=Chaudhary

 

 


[[alternative HTML version deleted]]

__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


RE: [R] Confidence Intervals from Bootstrap Replications

2004-12-17 Thread Berton Gunter

See , e.g. section 8.3 The two-sample problem of Efron and Tibshirani's AN
INTRODUCTION TO THE BOOTSTRAP. It makes it clear there why one just
bootstrap samples independently from the two separate samples.

The strata argument of boot() in the boot package allows one to do such
independent sampling and use the capabilities of that function.

-- Bert Gunter
Genentech Non-Clinical Statistics
South San Francisco, CA
 
The business of the statistician is to catalyze the scientific learning
process.  - George E. P. Box
 
 

 -Original Message-
 From: [EMAIL PROTECTED] 
 [mailto:[EMAIL PROTECTED] On Behalf Of 
 Mohammad A. Chaudhary
 Sent: Friday, December 17, 2004 2:39 PM
 To: [EMAIL PROTECTED]
 Subject: [R] Confidence Intervals from Bootstrap Replications 
 
 Hi All:
 
 I have to compute bootstrap confidence intervals, the statistic
 (incremental cost effectiveness ratio) is computed from two samples
 (intervention and control) of different sizes. All the bootstrap
 functions that I have seen use one dataset as argument. I may go ahead
 and get the desired number of bootstrap replications 
 separately. I would
 appreciate if you could point me to a source of a bootstrap 
 function (if
 available) that takes the B bootstrap replications and other 
 descriptive
 statistics and can get me the confidence intervals. Please 
 write me if I
 have not been clear in explaining my problem. Regards,
 
 Ashraf  
 
  
 
 ___
 
 M. Ashraf Chaudhary, Ph.D.
 
 Associate Scientist/Biostatistician
 
 Department of International Health
 
 Disease Prevention and Control Program
 
 Johns Hopkins University Bloomberg School of Public Health
 
 615 North Wolfe Street, Room W5506
 
 Baltimore MD 21205
 
  
 
 [EMAIL PROTECTED] 
 
 Phone: (410) 502-0741/Fax: (410) 502-6733
 
 http://faculty.jhsph.edu/?F=MohammadL=Chaudhary
 
  
 
  
 
 
   [[alternative HTML version deleted]]
 
 __
 [EMAIL PROTECTED] mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! 
 http://www.R-project.org/posting-guide.html


__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


RE: [R] confidence intervals

2004-09-09 Thread Rogers, James A [PGRD Groton]

Robert, 

I have this quick hack to obtain approximate Shewhart prediction intervals
for variance component models fit with lme (to nitpick slightly,
confidence intervals have the interpretation of containing parameters,
while prediction and tolerance intervals have the interpretation of
containing future observations or statistics). 

Back of the envelope documentation: the only argument that probably needs
explaining is the reps argument to shewhart(). If your model is, e.g.

fixed = Y ~ 1, random = ~ 1 | Batch

then specify reps = c(1, 1) if you want to predict a single future
observation from a single future batch, reps = c(1, 2) if you want to
predict the mean of two future observations from a single future batch, reps
= c(2, 2) if you want to predict the mean of 4 observations spread evenly
over 2 future batches, ...

Leave mult.check = 1, unless you want to do a Bonferroni correction. 

HTH, 

Jim Rogers


valStats2 - 
function (x, fixed, random, ...) 
{
mod - lme(fixed = fixed, data = x, random = random, ...)
mn - fixef(mod)
vc - VarCorr(mod)
err - Expecting only random intercept terms and a single fixed
intercept.\n
if (length(mn)  1 || ncol(vc)  2) 
stop(err)
rn - rownames(vc)
skip - regexpr(=, rn)  0
if (!any(skip)) 
vnms - attr(vc, title)
else vnms - grep(=, rn, value = TRUE)
vc - vc[!skip, ]
vnms - trim(sub(=.*, , vnms))
vnms - c(vnms, Residual)
vnms - paste(V, vnms, sep = .)
vars - as.numeric(vc[, Variance])
stats - c(mn, vars)
names(stats) - c(Intercept, vnms)
stats
}

shewhart - 
function (x, meancol = Intercept, varcols = grep(^V\\., names(x), value
= TRUE), reps = c(1, 1), alpha = 0.02, mult.check = 1) 
{
mn - x[[meancol]]
vr - as.matrix(x[varcols])
totvar - vr %*% (1/reps)
totsd - sqrt(totvar)
LL.mean - mn + qnorm(alpha/2/mult.check) * totsd
UL.mean - mn + qnorm(1 - alpha/2/mult.check) * totsd
out - data.frame(V.Total = totvar, LL.mean = LL.mean, UL.mean =
UL.mean)
out
}


### Example, where x is your data.frame:

foo - valStats2(x, fixed = Y ~ 1, random = ~ 1|Batch)
foo - as.data.frame(t(as.matrix(foo)))
data.frame(foo, shewhart(foo))



 -Original Message-
 From: [EMAIL PROTECTED]
 [mailto:[EMAIL PROTECTED] Behalf Of Spencer Graves
 Sent: Friday, September 03, 2004 11:58 AM
 To: [EMAIL PROTECTED]
 Cc: Robert Waters; [EMAIL PROTECTED]
 Subject: Re: [R] confidence intervals
 
 
 Hi, Robert: 
 
   While it may be difficult to program this in general 
 (as suggested 
 by it's position on Doug's To Do list), all the pieces should be 
 available to support a special script for your specific application.  
 What fixed and random model(s) interest you most? 
 
   hope this helps.  spencer graves
 
 Douglas Bates wrote:
 
  Robert Waters wrote:
 
  Dear R users;
 
  Im working with lme and Id like to have an idea of how
  can I get CI for the predictions made with the model.
  Im not a stats guy but, if Im not wrong, the CIs
  should be different if Im predicting a new data point
  or a new group. Ive been searching through the web and
  in help-lists with no luck. I know this topic had been
  asked before but without replies. Can anyone give an
  idea of where can I found information about this or
  how can I get it from R?
 
  Thanks for any hint
 
 
  That's not currently implemented in lme.  It's on the To 
 Do list but 
  it is not very close to the top.
 




LEGAL NOTICE\ Unless expressly stated otherwise, this messag...{{dropped}}

__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


RE: [R] confidence intervals

2004-09-08 Thread Warnes, Gregory R

You should be able to make small modifications to the ci.lme function
provided in the gregmisc/gmodels package.

-Greg

 -Original Message-
 From: [EMAIL PROTECTED]
 [mailto:[EMAIL PROTECTED] Behalf Of Spencer Graves
 Sent: Friday, September 03, 2004 11:58 AM
 To: [EMAIL PROTECTED]
 Cc: Robert Waters; [EMAIL PROTECTED]
 Subject: Re: [R] confidence intervals
 
 
 Hi, Robert: 
 
   While it may be difficult to program this in general 
 (as suggested 
 by it's position on Doug's To Do list), all the pieces should be 
 available to support a special script for your specific application.  
 What fixed and random model(s) interest you most? 
 
   hope this helps.  spencer graves
 
 Douglas Bates wrote:
 
  Robert Waters wrote:
 
  Dear R users;
 
  Im working with lme and Id like to have an idea of how
  can I get CI for the predictions made with the model.
  Im not a stats guy but, if Im not wrong, the CIs
  should be different if Im predicting a new data point
  or a new group. Ive been searching through the web and
  in help-lists with no luck. I know this topic had been
  asked before but without replies. Can anyone give an
  idea of where can I found information about this or
  how can I get it from R?
 
  Thanks for any hint
 
 
  That's not currently implemented in lme.  It's on the To 
 Do list but 
  it is not very close to the top.
 
  __
  [EMAIL PROTECTED] mailing list
  https://stat.ethz.ch/mailman/listinfo/r-help
  PLEASE do read the posting guide! 
  http://www.R-project.org/posting-guide.html
 
 
 -- 
 Spencer Graves, PhD, Senior Development Engineer
 O:  (408)938-4420;  mobile:  (408)655-4567
 
 __
 [EMAIL PROTECTED] mailing list
 https://stat.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! 
 http://www.R-project.org/posting-guide.html
 


LEGAL NOTICE\ Unless expressly stated otherwise, this messag...{{dropped}}

__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] confidence intervals

2004-09-03 Thread Douglas Bates
Robert Waters wrote:
Dear R users;
Im working with lme and Id like to have an idea of how
can I get CI for the predictions made with the model.
Im not a stats guy but, if Im not wrong, the CIs
should be different if Im predicting a new data point
or a new group. Ive been searching through the web and
in help-lists with no luck. I know this topic had been
asked before but without replies. Can anyone give an
idea of where can I found information about this or
how can I get it from R?
Thanks for any hint
That's not currently implemented in lme.  It's on the To Do list but 
it is not very close to the top.

__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] confidence intervals

2004-09-02 Thread Robert Waters
Dear R users;

Im working with lme and Id like to have an idea of how
can I get CI for the predictions made with the model.
Im not a stats guy but, if Im not wrong, the CIs
should be different if Im predicting a new data point
or a new group. Ive been searching through the web and
in help-lists with no luck. I know this topic had been
asked before but without replies. Can anyone give an
idea of where can I found information about this or
how can I get it from R?

Thanks for any hint

RW

__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] confidence intervals for linear combinations when using lme

2004-07-23 Thread Anna . H . Persson
Hi 

I really hope someone can help me.

I have just started to work with S-plus, and have not yet understood how it
really works. I am now trying to fit a mixed effects model with lme. My goal
is to compare four different groups, at several different time points, and I
therefore would like to create confidence intervals for linear combinations
of my estimated parameters (as I usually do with contrast or estimate or
lsmeans in SAS). I have now found out that the vcov-function in the
lme4-package could be of great help for me, but I do not know if I can use
lme4 when working with S-PLUS 6.1, or is it only for R. 

When I write 
library(MASS)
I find lme3, where there is a vcov-function, but if I understand it
correctly it works with lm and not with lme...

If I can use lme4 with s-plus 6.1, what shall I do in order to install it on
my computer?

I am grateful for all help I can get?

best regards
Anna Persson

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] confidence intervals for linear combinations when using lme

2004-07-23 Thread Deepayan Sarkar
On Friday 23 July 2004 05:36, [EMAIL PROTECTED] wrote:
 Hi

 I really hope someone can help me.

 I have just started to work with S-plus, and have not yet understood
 how it really works. I am now trying to fit a mixed effects model
 with lme. My goal is to compare four different groups, at several
 different time points, and I therefore would like to create
 confidence intervals for linear combinations of my estimated
 parameters (as I usually do with contrast or estimate or
 lsmeans in SAS). I have now found out that the vcov-function in the
 lme4-package could be of great help for me, but I do not know if I
 can use lme4 when working with S-PLUS 6.1, or is it only for R.

As far as I know, there are no plans (by the developers) to port lme4 to 
work with S-PLUS.

Deepayan

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] confidence intervals for linear combinations when using lme

2004-07-23 Thread Douglas Bates
[EMAIL PROTECTED] wrote:
I have just started to work with S-plus, and have not yet understood how it
really works. I am now trying to fit a mixed effects model with lme. My goal
is to compare four different groups, at several different time points, and I
therefore would like to create confidence intervals for linear combinations
of my estimated parameters (as I usually do with contrast or estimate or
lsmeans in SAS). I have now found out that the vcov-function in the
lme4-package could be of great help for me, but I do not know if I can use
lme4 when working with S-PLUS 6.1, or is it only for R. 
The lme4 package is only available for R at this time.
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Confidence intervals for predicted values in nls

2004-06-07 Thread joerg van den hoff
Cristina Silva wrote:
Dear all
I have tried to estimate the confidence intervals for predicted values of a
nonlinear model fitted with nls. The function predict gives the predicted
values and the lower and upper limits of the prediction, when the class of
the object is lm or glm. When the object is derived from nls, the function
predict (or predict.nls) gives only the predicted values. The se.fit and
interval aguments are just ignored.
Could anybody tell me how to estimate the confidence intervals for the
predicted values (not the model parameters), using an object of class nls?
Regards,
Cristina
--
Cristina Silva
IPIMAR - Departamento de Recursos Marinhos
Av. de Brasília, 1449-006 Lisboa
Portugal
Tel.: 351 21 3027096
Fax: 351 21 3015948
[EMAIL PROTECTED]
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
 

maybe this example helps:
==cut here===
#define a model formula (a and b are the parameters, f is x):
frml - k1 ~ f*(1-a*exp(-b/f))
#simulate some data:
a0 - .6
b0 - 1.2
f  - seq(0.01,4,length=20)
k1true- f*(1-a0*exp(-b0/f))
#include some noise
amp - .1
k1 - rnorm(k1true,k1true,amp*k1true)
#fit:
fifu - deriv(frml,c(a,b),function(a,b,x){})
rr-nls(k1~fifu(a,b,f),start=list(a=.5,b=1))
#the derivatives and variance/covariance matrix:
#(derivs could be extracted from fifu, too)
dk1.da - D(frml[[3]],'a')
dk1.db - D(frml[[3]],'b')
covar - vcov(rr)
#gaussian error propagation:
a - coef(rr)['a']
b - coef(rr)['b']
vark1 -
  eval(dk1.da)^2*covar[1,1]+
  eval(dk1.db)^2*covar[2,2]+
  2*eval(dk1.da)*eval(dk1.db)*covar[1,2]
errk1 - sqrt(vark1)
lower.bound - fitted(rr)-2*errk1  #two sigma !
upper.bound - fitted(rr)+2*errk1  #dito
plot(f,k1,pch=1)
ff - outer(c(1,1),f)
kk - outer(c(1,1),k1)*c(1-amp,1+amp)
matlines(ff,kk,lty=3,col=1)
matlines(f,cbind(k1true,fitted(rr),lower.bound,upper.bound),col=c(1,2,3,3),lty=c(1,1,2,2))
xylim - par('usr')
xpos - .1*(xylim[2]-xylim[1])+xylim[1]
ypos - xylim[4] - .1*(xylim[4]-xylim[3])
legend(xpos,ypos,
 c( 
'data',
'true',
'fit', 
'confidence'
  ), 
  pch=c(1,-1,-1,-1),
  lty=c(0,1,1,2),
  col=c(1,1,2,3)
)
==cut here===
if you put this in a file and source it a few times from within R you'll 
get an impression how often the fit deviates from the 'true' curve more 
than expected from
the shown confidence limits.

I believe this approach is 'nearly' valid as long as gaussian error 
probagation is valid (i.e. only to first order in covar and therefore 
for not too large std. errors, am I right?).
to my simple physicist's mind this should suffice to get 'reasonable' 
(probably, in strict sense, not completely correct?) confidence 
intervals for the fit/the prediction.
If somebody objects, please let me know!

joerg
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Confidence intervals for predicted values in nls

2004-06-04 Thread Alberto Murta
Parece que não tiveste grandes ajudas. Eu se fosse a ti fazia isso com um 
bootstrap bem planeado

Alberto


On Thursday 03 June 2004 21:22, Prof Brian Ripley wrote:
 On Thu, 3 Jun 2004, Cristina Silva wrote:
  I have tried to estimate the confidence intervals for predicted values of
  a nonlinear model fitted with nls. The function predict gives the
  predicted values and the lower and upper limits of the prediction, when
  the class of the object is lm or glm. When the object is derived from
  nls, the function predict (or predict.nls) gives only the predicted
  values. The se.fit and interval aguments are just ignored.

 Thre are no such arguments either to the generic function nls() nor its
 nls method.  Please do read the documentation!

  Could anybody tell me how to estimate the confidence intervals for the
  predicted values (not the model parameters), using an object of class
  nls?

 First you need to understand how to do this in theory: thereafter it is a
 programming task.  Hint: to find a confidence region for the parameters is
 not at all easy, as the examples in MASS (the book) and elsewhere show,
 and there is no guarantee that the confidence region for the prediction
 will be a single interval.

-- 
 Alberto G. Murta
Institute for Agriculture and Fisheries Research (INIAP-IPIMAR) 
Av. Brasilia, 1449-006 Lisboa, Portugal | Phone: +351 213027062
Fax:+351 213015948 | http://ipimar-iniap.ipimar.pt/pelagicos/

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Confidence intervals for predicted values in nls

2004-06-04 Thread Alberto Murta
Ops! 
I apologise for posting my last message to the R list by mistake.


-- 
 Alberto G. Murta
Institute for Agriculture and Fisheries Research (INIAP-IPIMAR) 
Av. Brasilia, 1449-006 Lisboa, Portugal | Phone: +351 213027062
Fax:+351 213015948 | http://ipimar-iniap.ipimar.pt/pelagicos/

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] Confidence intervals for predicted values in nls

2004-06-03 Thread Cristina Silva
Dear all

I have tried to estimate the confidence intervals for predicted values of a
nonlinear model fitted with nls. The function predict gives the predicted
values and the lower and upper limits of the prediction, when the class of
the object is lm or glm. When the object is derived from nls, the function
predict (or predict.nls) gives only the predicted values. The se.fit and
interval aguments are just ignored.

Could anybody tell me how to estimate the confidence intervals for the
predicted values (not the model parameters), using an object of class nls?

Regards,

Cristina

--
Cristina Silva
IPIMAR - Departamento de Recursos Marinhos
Av. de Brasília, 1449-006 Lisboa
Portugal
Tel.: 351 21 3027096
Fax: 351 21 3015948
[EMAIL PROTECTED]

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Confidence intervals for predicted values in nls

2004-06-03 Thread Prof Brian Ripley
On Thu, 3 Jun 2004, Cristina Silva wrote:

 I have tried to estimate the confidence intervals for predicted values of a
 nonlinear model fitted with nls. The function predict gives the predicted
 values and the lower and upper limits of the prediction, when the class of
 the object is lm or glm. When the object is derived from nls, the function
 predict (or predict.nls) gives only the predicted values. The se.fit and
 interval aguments are just ignored.

Thre are no such arguments either to the generic function nls() nor its 
nls method.  Please do read the documentation!

 Could anybody tell me how to estimate the confidence intervals for the
 predicted values (not the model parameters), using an object of class nls?

First you need to understand how to do this in theory: thereafter it is a 
programming task.  Hint: to find a confidence region for the parameters is 
not at all easy, as the examples in MASS (the book) and elsewhere show, 
and there is no guarantee that the confidence region for the prediction 
will be a single interval.

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Confidence intervals pointwise and family

2004-05-01 Thread Spencer Graves
?predict.lm in R 1.9.0 provides examples. 

hope this helps.  spencer graves
Linda portman wrote:
How can I add confidence intervals (pointwise and family) around the curves?
The curve is made by plot(x,y). Thanks!


-
[[alternative HTML version deleted]]
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
 

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Confidence Intervals for slopes

2004-04-01 Thread David Orme
Hi,

Thanks for the pointer - the 'lm(y~x:z)'  model does give the slopes 
directly and hence confint gives the confidence intervals.

The thing that puzzles me is that my dummy data explicitly sets the 
three levels of the factor to have different variances and yet the 
standard error is the same for all three parameter estimates in the 
summary.lm output - is this a common standard error of the 'x:z' term 
in the model? If you fit a separate regression to subsets of the data 
for each level in 'z' then the standard errors of the slope reflect 
these differences in variance. What I was trying to get was confidence 
limits from within a single model that also reflect the difference in 
certainty about the three slopes.

I realize that this is a failing of my understanding and more a stats 
question than an R question - if anyone can give me any advice or 
pointers, that would be great.

Thanks,
David
On 29 Mar 2004, at 20:04, BXC (Bendix Carstensen) wrote:

You may want:

lm( y ~ x:z )

This is the same model you fitted, but prametrized differently.
But please check that what you REALLY want is not
lm( y ~ z + x:z )

This is the model with different intercepts as well.

Bendix Carstensen
--
Bendix Carstensen
Senior Statistician
Steno Diabetes Center
Niels Steensens Vej 2
DK-2820 Gentofte
Denmark
tel: +45 44 43 87 38
mob: +45 30 75 87 38
fax: +45 44 43 07 06
[EMAIL PROTECTED]
www.biostat.ku.dk/~bxc
--




-Original Message-
From: [EMAIL PROTECTED]
[mailto:[EMAIL PROTECTED] On Behalf Of David Orme
Sent: Monday, March 29, 2004 4:44 PM
To: [EMAIL PROTECTED]
Subject: [R] Confidence Intervals for slopes
Hi,

I'm trying to get confidence intervals to slopes from a linear model
and I can't figure out how to get at them. As a cut 'n' paste example:
#
# dummy dataset - regression data for 3 treatments, each
treatment with
different (normal) variance
x - rep(1:10, length=30)
y - 10 - (rep(c(0.2,0.5,0.8), each=10)*x)+c(rnorm(10, sd=0.1),
rnorm(10, sd=0.6),rnorm(10, sd=1.1))
z - gl(3,10)
plot(y~x, pch=unclass(z))
# model as three slopes with common intercept
options(contrasts=c(contr.treatment,contr.poly))
model - lm(y~x+x:z)
# coefficient table in summary gives the intercept, first
slope and the
difference in slopes
summary(model)
# confint gives the confidence interval for the intercept and first
slope,
# and the CIs for the _differences_
confint(model)
#
What I'd like to report are the actual CI's for the slopes for the
second and third treatments, in the same way that confint returns the
parameter estimates for the first treatment. Can anyone point
me in the
right direction?
Thanks,
David
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo /r-help
PLEASE
do read the posting guide!
http://www.R-project.org/posting-guide.html

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] Confidence Intervals for slopes

2004-03-29 Thread David Orme
Hi,

I'm trying to get confidence intervals to slopes from a linear model 
and I can't figure out how to get at them. As a cut 'n' paste example:

#
# dummy dataset - regression data for 3 treatments, each treatment with 
different (normal) variance
x - rep(1:10, length=30)
y - 10 - (rep(c(0.2,0.5,0.8), each=10)*x)+c(rnorm(10, sd=0.1), 
rnorm(10, sd=0.6),rnorm(10, sd=1.1))
z - gl(3,10)
plot(y~x, pch=unclass(z))

# model as three slopes with common intercept
options(contrasts=c(contr.treatment,contr.poly))
model - lm(y~x+x:z)
# coefficient table in summary gives the intercept, first slope and the 
difference in slopes
summary(model)

# confint gives the confidence interval for the intercept and first 
slope,
# and the CIs for the _differences_
confint(model)
#

What I'd like to report are the actual CI's for the slopes for the 
second and third treatments, in the same way that confint returns the 
parameter estimates for the first treatment. Can anyone point me in the 
right direction?

Thanks,
David
__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Confidence Intervals for slopes

2004-03-29 Thread Andrew Robinson
David,

try the estimable() function in the gregmisc package.

Andrew

On Monday 29 March 2004 06:44, David Orme wrote:
 Hi,

 I'm trying to get confidence intervals to slopes from a linear model
 and I can't figure out how to get at them. As a cut 'n' paste example:

 #
 # dummy dataset - regression data for 3 treatments, each treatment with
 different (normal) variance
 x - rep(1:10, length=30)
 y - 10 - (rep(c(0.2,0.5,0.8), each=10)*x)+c(rnorm(10, sd=0.1),
 rnorm(10, sd=0.6),rnorm(10, sd=1.1))
 z - gl(3,10)
 plot(y~x, pch=unclass(z))

 # model as three slopes with common intercept
 options(contrasts=c(contr.treatment,contr.poly))
 model - lm(y~x+x:z)

 # coefficient table in summary gives the intercept, first slope and the
 difference in slopes
 summary(model)

 # confint gives the confidence interval for the intercept and first
 slope,
 # and the CIs for the _differences_
 confint(model)
 #

 What I'd like to report are the actual CI's for the slopes for the
 second and third treatments, in the same way that confint returns the
 parameter estimates for the first treatment. Can anyone point me in the
 right direction?

 Thanks,
 David

 __
 [EMAIL PROTECTED] mailing list
 https://www.stat.math.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide!
 http://www.R-project.org/posting-guide.html

-- 
Andrew Robinson  Ph: 208 885 7115
Department of Forest Resources   Fa: 208 885 6226
University of Idaho  E : [EMAIL PROTECTED]
PO Box 441133W : http://www.uidaho.edu/~andrewr
Moscow ID 83843  Or: http://www.biometrics.uidaho.edu
No statement above necessarily represents my employer's opinion.

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


RE: [R] Confidence Intervals for slopes

2004-03-29 Thread BXC (Bendix Carstensen)
You may want:

lm( y ~ x:z )

This is the same model you fitted, but prametrized differently.
But please check that what you REALLY want is not
  
lm( y ~ z + x:z )

This is the model with different intercepts as well.

Bendix Carstensen
--
Bendix Carstensen
Senior Statistician
Steno Diabetes Center
Niels Steensens Vej 2
DK-2820 Gentofte
Denmark
tel: +45 44 43 87 38
mob: +45 30 75 87 38
fax: +45 44 43 07 06
[EMAIL PROTECTED]
www.biostat.ku.dk/~bxc
--





 -Original Message-
 From: [EMAIL PROTECTED] 
 [mailto:[EMAIL PROTECTED] On Behalf Of David Orme
 Sent: Monday, March 29, 2004 4:44 PM
 To: [EMAIL PROTECTED]
 Subject: [R] Confidence Intervals for slopes
 
 
 Hi,
 
 I'm trying to get confidence intervals to slopes from a linear model 
 and I can't figure out how to get at them. As a cut 'n' paste example:
 
 #
 # dummy dataset - regression data for 3 treatments, each 
 treatment with 
 different (normal) variance
 x - rep(1:10, length=30)
 y - 10 - (rep(c(0.2,0.5,0.8), each=10)*x)+c(rnorm(10, sd=0.1), 
 rnorm(10, sd=0.6),rnorm(10, sd=1.1))
 z - gl(3,10)
 plot(y~x, pch=unclass(z))
 
 # model as three slopes with common intercept
 options(contrasts=c(contr.treatment,contr.poly))
 model - lm(y~x+x:z)
 
 # coefficient table in summary gives the intercept, first 
 slope and the 
 difference in slopes
 summary(model)
 
 # confint gives the confidence interval for the intercept and first 
 slope,
 # and the CIs for the _differences_
 confint(model)
 #
 
 What I'd like to report are the actual CI's for the slopes for the 
 second and third treatments, in the same way that confint returns the 
 parameter estimates for the first treatment. Can anyone point 
 me in the 
 right direction?
 
 Thanks,
 David
 
 __
 [EMAIL PROTECTED] mailing list 
 https://www.stat.math.ethz.ch/mailman/listinfo /r-help
 PLEASE 
 do read the posting guide! 
 http://www.R-project.org/posting-guide.html


__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] Confidence intervals for logistic regression

2004-02-19 Thread Michael Levine
Hi,

I found myself trying to figure out the type of confidence interval used for the 
coefficients of the logistic regression fit by using
glm(family=binomial)...

I suspect it is Wald confidence interval but am not sure...Does anybody know?
Also, if so, how can I ask for likelihood ratio and/or score-based confidence 
intervals?

Yours,
Michael
~
Michael Levine 
Assistant Professor
Department of Statistics
School of Arts and Sciences
Purdue University
Office: MATH 438
150 N. University Street
West Lafayette, IN 47907
phone  (765)496-7571
e-mail: [EMAIL PROTECTED]

[[alternative HTML version deleted]]

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] Confidence intervals for logistic regression

2004-02-19 Thread Prof Brian Ripley
On Thu, 19 Feb 2004, Michael Levine wrote:

 I found myself trying to figure out the type of confidence interval used
 for the coefficients of the logistic regression fit by using
 glm(family=binomial)...

AFAIK, no confidence interval is given by that call: it does not even 
calculate standard errors for the coefs (the summary method does that).
I could guess what you meant, but the answer depends on the guess 

 I suspect it is Wald confidence interval but am not sure...Does anybody
 know? Also, if so, how can I ask for likelihood ratio and/or score-based
 confidence intervals?

You can use confint() (whose glm method is in package MASS, which must be
attached) to give you profile-likelihood confidence intervals.  There
would be no point in trying to base confidence intervals on score tests,
as you would have to re-fit the model at all possible alternative values
and so the profile likelihood would be available.

-- 
Brian D. Ripley,  [EMAIL PROTECTED]
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel:  +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax:  +44 1865 272595

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


RE: [R] confidence-intervals in barchart

2004-02-10 Thread Marwan Khawaja
Try 'parplot2' -- 'gregmisc' package.
Marwan


---
Marwan Khawaja http://staff.aub.edu.lb/~mk36/
---


 -Original Message-
 From: [EMAIL PROTECTED]
 [mailto:[EMAIL PROTECTED] Behalf Of
 [EMAIL PROTECTED]
 Sent: Tuesday, February 10, 2004 4:01 AM
 To: [EMAIL PROTECTED]
 Subject: [R] confidence-intervals in barchart


 Hi R users,

 1)  How does one show confidence-intervals in a barchart and use rownames for
 labels on the y-axes?  I have looked at plotCI in gregmisc package . But
 it does not seem to produce something like a barchart.  The statistic, error,
 upper-bound, and lower-bound are in a dataframe.

 2) How to show CI in a barchart either using the statistic and, either (a)
 errors or (b) upper and lower bounds from a dataframe?

 An example will be very helpful.

   [[alternative HTML version deleted]]

 __
 [EMAIL PROTECTED] mailing list
 https://www.stat.math.ethz.ch/mailman/listinfo/r-help
 PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] confidence-intervals in barchart

2004-02-10 Thread Dominique Couturier
hello,
you can find a very detailed example of barplot with CI in the volume 
3/2 of R-news (october 2003).
hope this help.
dlc

Hi R users,

1)  How does one show confidence-intervals in a barchart and use 
rownames for
labels on the y-axes?  I have looked at plotCI in gregmisc 
package . But
it does not seem to produce something like a barchart.  The 
statistic, error,
upper-bound, and lower-bound are in a dataframe.

2) How to show CI in a barchart either using the statistic and, 
either (a)
errors or (b) upper and lower bounds from a dataframe?

An example will be very helpful.

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] confidence-intervals in barchart

2004-02-10 Thread Uwe Ligges
[EMAIL PROTECTED] wrote:

Hi R users,

1)  How does one show confidence-intervals in a barchart and use rownames for 
labels on the y-axes?  I have looked at plotCI in gregmisc package . But 
it does not seem to produce something like a barchart.  The statistic, error, 
upper-bound, and lower-bound are in a dataframe.

2) How to show CI in a barchart either using the statistic and, either (a) 
errors or (b) upper and lower bounds from a dataframe?

An example will be very helpful.

	[[alternative HTML version deleted]]

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


See the R Help Desk by Marc Schwartz in R News 3/2, 2003, pp. 2-6,
http://cran.r-project.org/doc/Rnews/Rnews_2003-2.pdf
Uwe Ligges

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] confidence-intervals in dotchart

2004-02-10 Thread TyagiAnupam
My earlier posting should have said dotchart, not barchart.

1)  How does one show confidence-intervals in a dotchart and use rownames 
for labels on the y-axes?  I have looked at plotCI in gregmisc package . 
But it does not seem to produce something like a dotchart.  The statistic, 
error, upper-bound, and lower-bound are in a dataframe.

2) How to show CI in a dotchart either using the statistic and, either (a) 
errors or (b) upper and lower bounds from a dataframe?

An example will be very helpful.


[[alternative HTML version deleted]]

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] confidence-intervals in dotchart

2004-02-10 Thread Frank E Harrell Jr
On Tue, 10 Feb 2004 04:30:46 EST
[EMAIL PROTECTED] wrote:

 My earlier posting should have said dotchart, not barchart.
 
 1)  How does one show confidence-intervals in a dotchart and use
 rownames for labels on the y-axes?  I have looked at plotCI in
 gregmisc package . But it does not seem to produce something like a
 dotchart.  The statistic, error, upper-bound, and lower-bound are in a
 dataframe.
 
 2) How to show CI in a dotchart either using the statistic and, either
 (a) errors or (b) upper and lower bounds from a dataframe?
 
 An example will be very helpful.

You might look at the example on p.149
http://cran.r-project.org/doc/contrib/Harrell-statcomp-notes.pdf and other
examples in that section.  These use the xYplot and Dotplot functions in
the Hmisc package, which extend lattice graphics to easily handle error
bars and bands.

---
Frank E Harrell Jr   Professor and Chair   School of Medicine
 Department of Biostatistics   Vanderbilt University

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


Re: [R] confidence-intervals in dotchart

2004-02-10 Thread TyagiAnupam
In a message dated 2/10/04 3:33:07 AM Pacific Standard Time, 
[EMAIL PROTECTED] writes:

 On Tue, 10 Feb 2004 04:30:46 EST
 [EMAIL PROTECTED] wrote:
 
 My earlier posting should have said dotchart, not barchart.
 
 1)  How does one show confidence-intervals in a dotchart and use
 rownames for labels on the y-axes?  I have looked at plotCI in
 gregmisc package . But it does not seem to produce something like a
 dotchart.  The statistic, error, upper-bound, and lower-bound are in a
 dataframe.
 
 2) How to show CI in a dotchart either using the statistic and, either
 (a) errors or (b) upper and lower bounds from a dataframe?
 
 An example will be very helpful.
 
 You might look at the example on p.149
 http://cran.r-project.org/doc/contrib/Harrell-statcomp-notes.pdf and other
 examples in that section.  These use the xYplot and Dotplot functions in
 the Hmisc package, which extend lattice graphics to easily handle error
 bars and bands.
 

This is wonderful! Thanks. How do I change the default colors for the session 
and the colors for each Dotplot() individually?  I can't do that from within 
the Dotplot() function.


[[alternative HTML version deleted]]

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] confidence intervals

2004-01-22 Thread Erin Hodgess
Hi Yun Fan!

Have you looked at predict.lm?

It can give you confidence and prediction intervals.

Hope this helps!

Sincerely,
Erin Hodgess
mailto: [EMAIL PROTECTED]

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


[R] Confidence intervals in ANOVA

2003-12-08 Thread Karl Knoblick

Hallo! 


I have the a model with 3 time points, 2 treatments 
and N subjects. I can calculate an ANOVA but I can not

calculate the CI of the interaction term (time and 
treatment), which I need for a closer look at the 
effect of the treatment to the 3 time points. I do NOT

want to use lme because I can not manage it to 
reproduce text book examples (see my posting [R] lme: 
reproducing example Karl Knoblick (Tue 02 Dec 2003 - 
21:34:54 EST)). 


Here some sample data: 


# Data 
# 35 subjects 
ID-factor(rep(1:35,each=3)) 
TREAT-factor(c(rep(A, 60), rep(B, 45))) 
TIME-factor(rep(1:3, 35)) 
Y-numeric(length=105) 
set.seed(1234) 
Y-rnorm(105) 
# want to see an effect:
Y[TREAT==A  TIME==2]-Y[TREAT==A  TIME==2] - 1
DF-data.frame(Y, ID, TREAT, TIME) 


# 2 possible designs: 
# Design 1 with random term 
DF.aov1-aov(Y ~ TIME*TREAT + Error(TREAT:ID), 
data=DF) 
summary(DF.aov1) 
# Design 2 without random term 
DF.aov2-aov(Y ~ TIME*TREAT, data=DF) 
summary(DF.aov2) 


I am also not sure about the design - I think design 1

is more appropriate. 


What I have tried is to calculate the CI of the 
coefficients: 
confint(DF.aov1[[2]]) 
confint(DF.aov1[[3]]) 


(or: 
confint(DF.aov2) 
) 


But how can I get the CI for a concrete difference for

example between the treatments at time point 2? 


I really hope, sombody can help! 


Karl

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help


[R] Confidence intervals plot

2003-06-14 Thread Ramon Martínez Coscollà
Hi all!!

I am trying to plot several confidence intervals in a unique plot. That is, for each 
x, I have a confidence interval for a parameter related to x and I would like to plot 
them in the same plot, in order to compare them. The plot would look like some 
parallel vertical lines, each one corresponding to a x value. Their extrem points 
would be the confidence interval limits.

I do not know if I am clear enough. Anyway, thank you in advance.

Ramon.
[[alternate HTML version deleted]]

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help


Re: [R] Confidence intervals plot

2003-06-14 Thread Duncan Murdoch
On Sat, 14 Jun 2003 20:44:09 +0200, you wrote:

Hi all!!

I am trying to plot several confidence intervals in a unique plot. That is, for each 
x, I have a confidence interval for a parameter related to x and I would like to plot 
them in the same plot, in order to compare them. The plot would look like some 
parallel vertical lines, each one corresponding to a x value. Their extrem points 
would be the confidence interval limits.

I do not know if I am clear enough. Anyway, thank you in advance.

Suppose the upper limits are in U, the lower limits in L, and the x
values in X.  Then 

# set up the axes etc.
plot(X, U, ylim=range(c(L,U)), type='n')

segments(X, L, X, U)

will do what you describe.  You should also look at arrows(), in case
you want points or crossbars on the ends of the segments.

Duncan Murdoch

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help


RE: [R] Confidence intervals plot

2003-06-14 Thread Marc Schwartz
-Original Message-
From: [EMAIL PROTECTED] 
[mailto:[EMAIL PROTECTED] On Behalf Of Ramon 
Martínez Coscollà
Sent: Saturday, June 14, 2003 1:44 PM
To: [EMAIL PROTECTED]
Subject: [R] Confidence intervals plot


Hi all!!

I am trying to plot several confidence intervals in a unique 
plot. That is, for each x, I have a confidence interval for a 
parameter related to x and I would like to plot them in the 
same plot, in order to compare them. The plot would look like 
some parallel vertical lines, each one corresponding to a x 
value. Their extrem points would be the confidence interval limits.

I do not know if I am clear enough. Anyway, thank you in advance.

Ramon.


You have several options depending upon whether you simply want
vertical CI lines above and below xy data points or if you might want
a barplot with CI's.

You could draw the xy data points using plot() and then draw the CI's
yourself using either segments() or arrows() in the base package or
see plotCI(), plotmeans() and barplot2() in the 'gregmisc' package on
CRAN. plotmeans() is a wrapper function that can call plotCI().

HTH,

Marc Schwartz

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help


Re: [R] Confidence intervals plot

2003-06-14 Thread kjetil brinchmann halvorsen
On 14 Jun 2003 at 20:44, Ramon Martínez Coscollà wrote:

Hola!

First set up a plot with

plot( c(x.lower, x.upper) , c(y.lower, y.upper), type=n )

and then add each line using segments

Kjetil Halvorsen

 Hi all!!
 
 I am trying to plot several confidence intervals in a unique plot. That is, for each 
 x, I have a confidence interval for a parameter related to x and I would like to 
 plot them in the same plot, in order to compare them. The plot would look like some 
 parallel vertical lines, each one 
corresponding to a x value. Their extrem points would be the confidence interval 
limits.
 
 I do not know if I am clear enough. Anyway, thank you in advance.
 
 Ramon.
   [[alternate HTML version deleted]]
 
 __
 [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