[R] Confidence interval for relative risk

2006-11-13 Thread David Duffy
Michael Dewey [EMAIL PROTECTED] wrote
 Subject: [R] Confidence interval for relative risk

 The concrete problem is that I am refereeing
 a paper where a confidence interval is
 presented for the risk ratio and I do not find
 it credible. I show below my attempts to
 do this in R. The example is slightly changed
 from the authors'.

 I can obtain a confidence interval for
 the odds ratio from fisher.test of
 course

 === fisher.test example ===

  outcome - matrix(c(500, 0, 500, 8), ncol = 2, byrow = TRUE)
  fisher.test(outcome)

  Fisher's Exact Test for Count Data

 data:  outcome
 p-value = 0.00761
 alternative hypothesis: true odds ratio is not equal to 1
 95 percent confidence interval:
   1.694792  Inf
 sample estimates:
 odds ratio
 Inf

 === end example ===


Since

RR   = p1 +1   p1=a/(a+b), p2=c/(c+d)
 --
   p2 1
 -   +   ---
  1-p2OR

you can backsolve for the RR CI.

 x - function(p1, p2, oddsr) p1 + 1/(p2/(1-p2) + 1/oddsr)
 1/x(8/508, 0, 1/1.694792)
[1] 1.650734

Another approach is the Cornfield method -- which is a profile likelihood based 
CI 
using the Gibbs Chi-square (LRTS), but IIRC originally used the Pearson 
chi-square.
I don't think there are R implementations.

David Duffy.

-- 
| David Duffy (MBBS PhD) ,-_|\
| email: [EMAIL PROTECTED]  ph: INT+61+7+3362-0217 fax: -0101  / *
| Epidemiology Unit, Queensland Institute of Medical Research   \_,-._/
| 300 Herston Rd, Brisbane, Queensland 4029, Australia  GPG 4D0B994A v

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

2006-11-12 Thread Michael Dewey
At 15:54 10/11/2006, Viechtbauer Wolfgang (STAT) wrote:

Thanks for the suggestion Wolfgang, but whatever the original authors 
did that is not it.

Hello,

A common way to calculate the CI for a relative risk is the 
following. Given a 2x2 table of the form:

a b
c d

the log of the risk ratio is given by:

lrr = log[ a/(a+b) / ( c/(c+d) ) ]

which is asymptotically normal with variance:

vlrr = 1/a - 1/(a+b) + 1/c - 1/(c+d).

So an approximate 95% CI for the risk ratio is given by:

exp[ lrr - 1.96*sqrt(vlrr) ], exp[ lrr + 1.96*sqrt(vlrr) ].

A common convention is to add 1/2 to each cell when there are zeros.

So, for the table:

Col 1 Col 2
Row 1 8   500
Row 2 0   500

lrr = log[ 8.5/509 / ( 0.5/501 ) ] = 2.817
vllr = 1/8.5 - 1/509 + 1/0.5 - 1/501 = 2.1137

exp[ 2.817-1.96*sqrt(2.1137) ] = .97
exp[ 2.817+1.96*sqrt(2.1137) ] = 289.04

Maybe that is what the authors did.

Best,

--
Wolfgang Viechtbauer
  Department of Methodology and Statistics
  University of Maastricht, The Netherlands
  http://www.wvbauer.com/


  -Original Message-
  From: [EMAIL PROTECTED] [mailto:r-help-
  [EMAIL PROTECTED] On Behalf Of Michael Dewey
  Sent: Friday, November 10, 2006 15:43
  To: r-help@stat.math.ethz.ch
  Subject: [R] Confidence interval for relative risk
 
  The concrete problem is that I am refereeing
  a paper where a confidence interval is
  presented for the risk ratio and I do not find
  it credible. I show below my attempts to
  do this in R. The example is slightly changed
  from the authors'.
 
  I can obtain a confidence interval for
  the odds ratio from fisher.test of
  course
 
  === fisher.test example ===
 
outcome - matrix(c(500, 0, 500, 8), ncol = 2, byrow = TRUE)
fisher.test(outcome)
 
   Fisher's Exact Test for Count Data
 
  data:  outcome
  p-value = 0.00761
  alternative hypothesis: true odds ratio is not equal to 1
  95 percent confidence interval:
1.694792  Inf
  sample estimates:
  odds ratio
  Inf
 
  === end example ===
 
  but in epidemiology authors often
  prefer to present risk ratios.
 
  Using the facility on CRAN to search
  the site I find packages epitools and Epi
  which both offer confidence intervals
  for the risk ratio
 
  === Epi example ===
 
library(Epi)
twoby2(outcome[c(2,1),c(2,1)])
  2 by 2 table analysis:
  --
  Outcome   : Col 1
  Comparing : Row 1 vs. Row 2
 
 Col 1 Col 2P(Col 1) 95% conf. interval
  Row 1 8   500  0.01570.0079   0.0312
  Row 2 0   500  0.0.  NaN
 
  95% conf. interval
Relative Risk:Inf   NaN  Inf
Sample Odds Ratio:Inf   NaN  Inf
  Conditional MLE Odds Ratio:Inf1.6948  Inf
   Probability difference: 0.01570.0027   0.0337
 
Exact P-value: 0.0076
   Asymptotic P-value: NaN
  --
 
  === end example ===
 
  So Epi gives me a lower limit of NaN but the same confidence
  interval and p-value as fisher.test
 
  === epitools example ===
 
library(epitools)
riskratio(outcome)
  $data
 Outcome
  Predictor  Disease1 Disease2 Total
 Exposed1  5000   500
 Exposed2  5008   508
 Total10008  1008
 
  $measure
 risk ratio with 95% C.I.
  Predictor  estimate lower upper
 Exposed11NANA
 Exposed2  Inf   NaN   Inf
 
  $p.value
 two-sided
  Predictor  midp.exact fisher.exact  chi.square
 Exposed1 NA   NA  NA
 Exposed2 0.00404821  0.007610478 0.004843385
 
  $correction
  [1] FALSE
 
  attr(,method)
  [1] Unconditional MLE  normal approximation (Wald) CI
  Warning message:
  Chi-squared approximation may be incorrect in: chisq.test(xx, correct =
  correction)
 
  === end example ===
 
  And epitools also gives a lower limit
  of NaN.
 
  === end all examples ===
 
  I would prefer not to have to tell the authors of the
  paper I am refereeing that
  I think they are wrong unless I can help them with what they
  should have done.
 
  Is there another package I should have tried?
 
  Is there some other way of doing this?
 
  Am I doing something fundamentally wrong-headed?
 
 
 
  Michael Dewey
  http://www.aghmed.fsnet.co.uk

Michael Dewey
http://www.aghmed.fsnet.co.uk

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

2006-11-12 Thread Peter Dalgaard
Michael Dewey [EMAIL PROTECTED] writes:

 At 15:54 10/11/2006, Viechtbauer Wolfgang (STAT) wrote:
 
 Thanks for the suggestion Wolfgang, but whatever the original authors 
 did that is not it.

Did you ever say what result they got?

-p
 
 Hello,
 
 A common way to calculate the CI for a relative risk is the 
 following. Given a 2x2 table of the form:
 
 a b
 c d
 
 the log of the risk ratio is given by:
 
 lrr = log[ a/(a+b) / ( c/(c+d) ) ]
 
 which is asymptotically normal with variance:
 
 vlrr = 1/a - 1/(a+b) + 1/c - 1/(c+d).
 
 So an approximate 95% CI for the risk ratio is given by:
 
 exp[ lrr - 1.96*sqrt(vlrr) ], exp[ lrr + 1.96*sqrt(vlrr) ].
 
 A common convention is to add 1/2 to each cell when there are zeros.
 
 So, for the table:
 
 Col 1 Col 2
 Row 1 8   500
 Row 2 0   500
 
 lrr = log[ 8.5/509 / ( 0.5/501 ) ] = 2.817
 vllr = 1/8.5 - 1/509 + 1/0.5 - 1/501 = 2.1137
 
 exp[ 2.817-1.96*sqrt(2.1137) ] = .97
 exp[ 2.817+1.96*sqrt(2.1137) ] = 289.04
 
 Maybe that is what the authors did.
 
 Best,
 
 --
 Wolfgang Viechtbauer
   Department of Methodology and Statistics
   University of Maastricht, The Netherlands
   http://www.wvbauer.com/
 
 
   -Original Message-
   From: [EMAIL PROTECTED] [mailto:r-help-
   [EMAIL PROTECTED] On Behalf Of Michael Dewey
   Sent: Friday, November 10, 2006 15:43
   To: r-help@stat.math.ethz.ch
   Subject: [R] Confidence interval for relative risk
  
   The concrete problem is that I am refereeing
   a paper where a confidence interval is
   presented for the risk ratio and I do not find
   it credible. I show below my attempts to
   do this in R. The example is slightly changed
   from the authors'.
  
   I can obtain a confidence interval for
   the odds ratio from fisher.test of
   course
  
   === fisher.test example ===
  
 outcome - matrix(c(500, 0, 500, 8), ncol = 2, byrow = TRUE)
 fisher.test(outcome)
  
Fisher's Exact Test for Count Data
  
   data:  outcome
   p-value = 0.00761
   alternative hypothesis: true odds ratio is not equal to 1
   95 percent confidence interval:
 1.694792  Inf
   sample estimates:
   odds ratio
   Inf
  
   === end example ===
  
   but in epidemiology authors often
   prefer to present risk ratios.
  
   Using the facility on CRAN to search
   the site I find packages epitools and Epi
   which both offer confidence intervals
   for the risk ratio
  
   === Epi example ===
  
 library(Epi)
 twoby2(outcome[c(2,1),c(2,1)])
   2 by 2 table analysis:
   --
   Outcome   : Col 1
   Comparing : Row 1 vs. Row 2
  
  Col 1 Col 2P(Col 1) 95% conf. interval
   Row 1 8   500  0.01570.0079   0.0312
   Row 2 0   500  0.0.  NaN
  
   95% conf. interval
 Relative Risk:Inf   NaN  Inf
 Sample Odds Ratio:Inf   NaN  Inf
   Conditional MLE Odds Ratio:Inf1.6948  Inf
Probability difference: 0.01570.0027   0.0337
  
 Exact P-value: 0.0076
Asymptotic P-value: NaN
   --
  
   === end example ===
  
   So Epi gives me a lower limit of NaN but the same confidence
   interval and p-value as fisher.test
  
   === epitools example ===
  
 library(epitools)
 riskratio(outcome)
   $data
  Outcome
   Predictor  Disease1 Disease2 Total
  Exposed1  5000   500
  Exposed2  5008   508
  Total10008  1008
  
   $measure
  risk ratio with 95% C.I.
   Predictor  estimate lower upper
  Exposed11NANA
  Exposed2  Inf   NaN   Inf
  
   $p.value
  two-sided
   Predictor  midp.exact fisher.exact  chi.square
  Exposed1 NA   NA  NA
  Exposed2 0.00404821  0.007610478 0.004843385
  
   $correction
   [1] FALSE
  
   attr(,method)
   [1] Unconditional MLE  normal approximation (Wald) CI
   Warning message:
   Chi-squared approximation may be incorrect in: chisq.test(xx, correct =
   correction)
  
   === end example ===
  
   And epitools also gives a lower limit
   of NaN.
  
   === end all examples ===
  
   I would prefer not to have to tell the authors of the
   paper I am refereeing that
   I think they are wrong unless I can help them with what they
   should have done.
  
   Is there another package I should have tried?
  
   Is there some other way of doing this?
  
   Am I doing something fundamentally wrong-headed?
  
  
  
   Michael Dewey
   http://www.aghmed.fsnet.co.uk
 
 Michael Dewey
 http://www.aghmed.fsnet.co.uk
 
 __
 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

Re: [R] Confidence interval for relative risk

2006-11-12 Thread Michael Dewey
At 12:35 12/11/2006, Peter Dalgaard wrote:
Michael Dewey [EMAIL PROTECTED] writes:

  At 15:54 10/11/2006, Viechtbauer Wolfgang (STAT) wrote:
 
  Thanks for the suggestion Wolfgang, but whatever the original authors
  did that is not it.

Did you ever say what result they got?

 -p


No, because I did not want to use the original numbers in the 
request. So as the snippet below indicates I changed the numbers. If 
I apply Wolfgang's suggestion (which I had already thought of but 
discarded) I get about 13 for the real example where the authors quote about 5.

My question still remains though as to how I can get a confidence 
interval for the risk ratio without adding a constant to each cell.

it credible. I show below my attempts to
do this in R. The example is slightly changed
from the authors'.

Michael Dewey
http://www.aghmed.fsnet.co.uk

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

2006-11-12 Thread Spencer Graves
  When I have refereed manuscripts for publication and have been 
unable to get their answers, I have told the authors they need more 
explanation of their methodology -- while summarizing what I tried in a 
few lines.  I've even told some that it would make it vastly easier for 
the referees and increase the potential readership for their article if 
they make R code available -- at least downloadable somewhere and 
preferably in a contributed R package, where it could attract interest 
in the article from audiences who would not likely find it any other 
way.  I've done this for articles that were NOT written to specifically 
describe statistical software. 

  I have not followed all of this thread.  However, one suggestion I 
saw looked like it used the delta method, if I understood correctly 
from skimming without studying the details carefully.  Have you also 
considered 2*log(likelihood ratio) being approximately chi-square? 

  Just my 2e-9 Euros (or 2e-7 Yen or Yuan). 
   Spencer Graves

Michael Dewey wrote:
 At 12:35 12/11/2006, Peter Dalgaard wrote:
   
 Michael Dewey [EMAIL PROTECTED] writes:

 
 At 15:54 10/11/2006, Viechtbauer Wolfgang (STAT) wrote:

 Thanks for the suggestion Wolfgang, but whatever the original authors
 did that is not it.
   
 Did you ever say what result they got?

 -p

 

 No, because I did not want to use the original numbers in the 
 request. So as the snippet below indicates I changed the numbers. If 
 I apply Wolfgang's suggestion (which I had already thought of but 
 discarded) I get about 13 for the real example where the authors quote about 
 5.

 My question still remains though as to how I can get a confidence 
 interval for the risk ratio without adding a constant to each cell.

   
 it credible. I show below my attempts to
 do this in R. The example is slightly changed
 from the authors'.
   

 Michael Dewey
 http://www.aghmed.fsnet.co.uk

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

2006-11-10 Thread Michael Dewey
The concrete problem is that I am refereeing
a paper where a confidence interval is
presented for the risk ratio and I do not find
it credible. I show below my attempts to
do this in R. The example is slightly changed
from the authors'.

I can obtain a confidence interval for
the odds ratio from fisher.test of
course

=== fisher.test example ===

  outcome - matrix(c(500, 0, 500, 8), ncol = 2, byrow = TRUE)
  fisher.test(outcome)

 Fisher's Exact Test for Count Data

data:  outcome
p-value = 0.00761
alternative hypothesis: true odds ratio is not equal to 1
95 percent confidence interval:
  1.694792  Inf
sample estimates:
odds ratio
Inf

=== end example ===

but in epidemiology authors often
prefer to present risk ratios.

Using the facility on CRAN to search
the site I find packages epitools and Epi
which both offer confidence intervals
for the risk ratio

=== Epi example ===

  library(Epi)
  twoby2(outcome[c(2,1),c(2,1)])
2 by 2 table analysis:
--
Outcome   : Col 1
Comparing : Row 1 vs. Row 2

   Col 1 Col 2P(Col 1) 95% conf. interval
Row 1 8   500  0.01570.0079   0.0312
Row 2 0   500  0.0.  NaN

95% conf. interval
  Relative Risk:Inf   NaN  Inf
  Sample Odds Ratio:Inf   NaN  Inf
Conditional MLE Odds Ratio:Inf1.6948  Inf
 Probability difference: 0.01570.0027   0.0337

  Exact P-value: 0.0076
 Asymptotic P-value: NaN
--

=== end example ===

So Epi gives me a lower limit of NaN but the same confidence
interval and p-value as fisher.test

=== epitools example ===

  library(epitools)
  riskratio(outcome)
$data
   Outcome
Predictor  Disease1 Disease2 Total
   Exposed1  5000   500
   Exposed2  5008   508
   Total10008  1008

$measure
   risk ratio with 95% C.I.
Predictor  estimate lower upper
   Exposed11NANA
   Exposed2  Inf   NaN   Inf

$p.value
   two-sided
Predictor  midp.exact fisher.exact  chi.square
   Exposed1 NA   NA  NA
   Exposed2 0.00404821  0.007610478 0.004843385

$correction
[1] FALSE

attr(,method)
[1] Unconditional MLE  normal approximation (Wald) CI
Warning message:
Chi-squared approximation may be incorrect in: chisq.test(xx, correct =
correction)

=== end example ===

And epitools also gives a lower limit
of NaN.

=== end all examples ===

I would prefer not to have to tell the authors of the
paper I am refereeing that
I think they are wrong unless I can help them with what they
should have done.

Is there another package I should have tried?

Is there some other way of doing this?

Am I doing something fundamentally wrong-headed?



Michael Dewey
http://www.aghmed.fsnet.co.uk

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