#| 

  Cond exponential in Racket/Gamble 

  From PSI test/condExponential.psi
  """
  https://www.quora.com/If-X-and-Y-are-independent-exponential-random-variables-with-parameter-1-What-is-E-X+Y-2-X-
1-Y-1
  """
  $ psi --expectation condExponential.psi
  E(r_) = 18

  Mathematica:
  """
  Expectation((x + y)^2 | x > 1 && y > 1, (x ~ ExponentialDistribution(1), y ~ ExponentialDistribution(1)))
  -> 18
  """

var : x
mean: 1.9913873540132898
Min: 1.0000709398886531 Mean: 2.0035821043743915 Max: 11.387193288268117 Variance: 1.0138029278118805 Stddev: 1.0068778117586465
Credible interval (0.94): 1.0000709398886531..3.826930413474771
Percentiles:
(0.01 1.009995222874919)
(0.025 1.0235203751058912)
(0.1 1.1057226425112068)
(0.05 1.0481741983929291)
(0.25 1.2907512416507527)
(0.5 1.698870568732206)
(0.75 2.3895950867165445)
(0.84 2.8418239670792245)
(0.9 3.3108688683702114)
(0.95 4.000655883758724)
(0.975 4.6808414500531015)
(0.99 5.574298270391136)
(0.999 8.12873197548311)
Histogram:
1     : 1   
1.236 : 2101
1.472 : 1647
1.708 : 1306
1.944 : 1092
2.18  : 779 
2.416 : 645 
2.653 : 508 
2.889 : 401 
3.125 : 326 
3.361 : 245 
3.597 : 203 
3.833 : 153 
4.069 : 127 
4.305 : 92  
4.541 : 80  
4.777 : 62  
5.013 : 45  
5.249 : 36  
5.485 : 39  
5.721 : 25  
5.958 : 20  
6.194 : 9   
6.43  : 11  
6.666 : 13  
6.902 : 6   
7.138 : 3   
7.374 : 3   
7.61  : 5   
7.846 : 5   
8.082 : 1   
8.318 : 2   
8.554 : 3   
8.79  : 1   
9.026 : 0   
9.263 : 0   
9.499 : 1   
9.735 : 2   
9.971 : 1   
10.207: 0   
10.443: 0   
10.679: 0   
10.915: 0   
11.151: 0   

var : y
mean: 1.9927718095087743
Min: 1.0000108904741494 Mean: 2.003747496481179 Max: 10.23154426542632 Variance: 1.039756613441866 Stddev: 1.019684565658354
Credible interval (0.94): 1.0000108904741494..3.888236404266661
Percentiles:
(0.01 1.0098929616820953)
(0.025 1.02446603898925)
(0.1 1.1025718798084785)
(0.05 1.0519745537771819)
(0.25 1.281822585708286)
(0.5 1.6906569364599207)
(0.75 2.37858797539471)
(0.84 2.8452504083928503)
(0.9 3.326564082688571)
(0.95 4.066721401624673)
(0.975 4.735134914724011)
(0.99 5.735280256028153)
(0.999 8.253042601448374)
Histogram:
1     : 1   
1.21  : 1925
1.42  : 1532
1.629 : 1220
1.839 : 1038
2.049 : 809 
2.259 : 660 
2.469 : 546 
2.678 : 403 
2.888 : 329 
3.098 : 274 
3.308 : 241 
3.518 : 182 
3.728 : 165 
3.937 : 107 
4.147 : 105 
4.357 : 91  
4.567 : 71  
4.777 : 57  
4.986 : 43  
5.196 : 37  
5.406 : 30  
5.616 : 18  
5.826 : 28  
6.035 : 16  
6.245 : 13  
6.455 : 13  
6.665 : 7   
6.875 : 6   
7.084 : 3   
7.294 : 8   
7.504 : 4   
7.714 : 2   
7.924 : 3   
8.133 : 1   
8.343 : 2   
8.553 : 5   
8.763 : 0   
8.973 : 1   
9.183 : 1   
9.392 : 1   
9.602 : 1   
9.812 : 0   
10.022: 0   

var : e
mean: 17.853097746257713
Min: 4.023641724534582 Mean: 18.073831233876092 Max: 188.28294134783667 Variance: 216.6160114470597 Stddev: 14.71788067104295
Credible interval (0.94): 4.0635878252905915..42.98469422707255
Percentiles:
(0.01 4.628827774117357)
(0.025 5.046769407468217)
(0.1 6.462747020826413)
(0.05 5.581213537475119)
(0.25 8.812158914305126)
(0.5 13.637104712717898)
(0.75 22.05360473892758)
(0.84 28.082439025436905)
(0.9 34.83942252145432)
(0.95 45.175895595188685)
(0.975 57.07787354958004)
(0.99 76.42455229619361)
(0.999 129.41640644495376)
Histogram:
4.024  : 1   
8.211  : 2126
12.399 : 2307
16.587 : 1705
20.774 : 1115
24.962 : 742 
29.15  : 504 
33.338 : 378 
37.525 : 281 
41.713 : 184 
45.901 : 186 
50.088 : 109 
54.276 : 78  
58.464 : 53  
62.652 : 49  
66.839 : 36  
71.027 : 21  
75.215 : 21  
79.402 : 18  
83.59  : 16  
87.778 : 12  
91.966 : 8   
96.153 : 8   
100.341: 3   
104.529: 9   
108.716: 4   
112.904: 5   
117.092: 3   
121.28 : 2   
125.467: 1   
129.655: 5   
133.843: 1   
138.03 : 2   
142.218: 1   
146.406: 1   
150.594: 0   
154.781: 0   
158.969: 0   
163.157: 2   
167.344: 0   
171.532: 1   
175.72 : 0   
179.908: 0   
184.095: 1   



  This program was created by Hakan Kjellerstrand, hakank@gmail.com
  See also my Racket page: http://www.hakank.org/racket/

|#

#lang gamble

; (require gamble/viz)
(require racket)
(require "gamble_utils.rkt")
; (require "gamble_distributions.rkt")

(define (model1)
  (; enumerate
   ; rejection-sampler
   importance-sampler
   ; mh-sampler

   (define x (exponential 1))
   (define y (exponential 1))

   (observe/fail (and (> x 1) (> y 1)))
    
   (define e (expt (+ x y) 2))

   (list x
         y
         e
         )
   
   )
)

(show-marginals (model1)
                (list  "x"
                       "y"
                       "e"
                       )
                #:num-samples 10000
                #:truncate-output 5
                #:skip-marginals? #t
                #:show-stats? #t
                #:credible-interval 0.94
                #:show-histogram? #t
                #:show-percentiles? #t
                )