*!  Version 1.0.1 <JRG; 10Sep1999>   (STB-52 sg119)
program define propcii
*!  CIs for a proportion, immediate
*!  Syntax:  . prcii <n> <x> [, Level(0.xx)
*!                       ALL ARcsin EWald EXact WAld WIlson ]


   version 6.0
   if "`1'" == "?" {
      which propcii
      exit
   }

   gettoken Nn 0 : 0, parse(" ,")
   gettoken X  0 : 0, parse(" ,")
   confirm integer number `Nn'
   confirm integer number `X'
   if !( (`X'*(`Nn'-`X')>=0) * (`Nn'>0) ) {
      di in red "Invalid arguments: `Nn' `X'"
      error 499
   }
   syntax [, Level(real $S_level) ALL ARcsin EWald EXact WAld WIlson]

   if "`all'" != "" { local A "exact wald ewald wilson arcsin" }
   else {
      local A "`exact' `wald' `ewald' `wilson' `arcsin'"
      if "`A'" == "    " { local A "exact" }
   }
   if `level' < 1 { local level = 100*`level' }
   if `X' < 0.5*(`Nn'+1) {
      local X = `Nn' - `X'
      local UM "1-"
   }

   global Prop_Lub = ((100-`level')/200)^(1/`Nn')
   global Prop_Uub = ((100+`level')/200)^(1/`Nn')
   local p_hat = `X' / `Nn'
   x -1 `Nn' `UM'`p_hat' `level'
   global Prop___U 1
   global Prop___L = $Prop_Lub
   local z = invnorm((100+`level')/200)

   tokenize `A'
   while "`1'" != "" {
      if "`1'" == "exact" {
         global S_1 "Exact (Clopper-Pearson)"
         if `X' < `Nn' {
            global Prop___L = invbinomial(`Nn',`X',(100-`level')/200)
            global Prop___U = cond(`X' < `Nn'-1, /*
               */  invbinomial(`Nn',`X',(100+`level')/200),$Prop_Uub)
         }
         x $Prop___L $Prop___U S_1 `UM'
      }
      else if "`1'" == "wald" {
         global S_1 "   Wald (normal theory)"
         if `X' < `Nn' { x0 `X' `Nn' `z' 0 }
         x $Prop___L $Prop___U S_1 `UM'
      }
      else if "`1'" == "ewald" {
         global S_1 "  EWald (Enhanced Wald)"
         if `X' < `Nn' {
            local L = 2.64*( `p_hat'*(1-`p_hat') )^0.2
            x0 `X' `Nn' `z' `L'
         }
         x $Prop___L $Prop___U S_1 `UM'
      }
      else if "`1'" == "wilson" {
         global S_1 " Wilson (Agresti-Coull)"
         if `X' < `Nn' { x0 `X' `Nn' `z' 2 }
         x $Prop___L $Prop___U S_1 `UM'
      }
      else if "`1'" == "arcsin" {
         global S_1 " Arcsin transform-based"
         if `X' < `Nn' {
            local L = 0.25*`X'/`Nn'
            local p = asin(sqrt((`X'+`L')/(`Nn'+2*`L')))
            local SE = 0.5/sqrt(`Nn')
            global Prop___L = sin( max(0, `p' - `z'*`SE') )^2
            global Prop___U = cond(`X'<`Nn'-1,    /*
               */    sin( min(_pi/2, `p'+`z'*`SE') )^2, $Prop_Uub)
         }
         x $Prop___L $Prop___U S_1 `UM'
      }
      mac shift
   }
   mac drop Prop_*
end


program define x, rclass
   if `1' == -1 {
      return clear
      return scalar N = `2'
      return scalar p_hat = `3'
      return scalar level = `4'
      exit
   }

   if r(p_hat)*(1-r(p_hat))  {
      local 1 = max(1 - $Prop_Uub, min(`1', $Prop_Lub))
      local 2 = max(1 - $Prop_Lub, min(`2', $Prop_Uub))
   }
   if "`4'" != "" {
      local 4 = `1'
      local 1 = 1 - `2'
      local 2 = 1 - `4'
   }
   noi di in gr "${`3'}", r(level) "% CI: [" in ye %6.4f `1'   /*
         */    in gr ", " in ye  %6.4f `2' in gr "]"

   local a : word 1 of ${`3'}
   local a = substr("`a'", 1, 5)
   return scalar ub_`a' = `2'
   return scalar lb_`a' = `1'
   return add
end


program define x0
   local p = (`1'+`4')/(`2'+2*`4')
   local SE = sqrt(`p'*(1-`p')/(`2'+2*`4'))
   global Prop___L = `p' - `3'*`SE'
   global Prop___U = `p' + `3'*`SE'
end