src/HOL/Probability/Distributions.thy
changeset 57235 b0b9a10e4bf4
parent 56996 891e992e510f
child 57252 19b7ace1c5da
     1.1 --- a/src/HOL/Probability/Distributions.thy	Wed Jun 11 13:39:38 2014 +0200
     1.2 +++ b/src/HOL/Probability/Distributions.thy	Thu Jun 12 15:47:36 2014 +0200
     1.3 @@ -1,14 +1,313 @@
     1.4 +(*  Title:      HOL/Probability/Distributions.thy
     1.5 +    Author:     Sudeep Kanav, TU München
     1.6 +    Author:     Johannes Hölzl, TU München *)
     1.7 +
     1.8 +header {* Properties of Various Distributions *}
     1.9 +
    1.10  theory Distributions
    1.11 -  imports Probability_Measure
    1.12 +  imports Probability_Measure Convolution
    1.13  begin
    1.14  
    1.15 +lemma measure_lebesgue_Icc: "measure lebesgue {a .. b} = (if a \<le> b then b - a else 0)"
    1.16 +  by (auto simp: measure_def)
    1.17 +
    1.18 +lemma integral_power:
    1.19 +  "a \<le> b \<Longrightarrow> (\<integral>x. x^k * indicator {a..b} x \<partial>lborel) = (b^Suc k - a^Suc k) / Suc k"
    1.20 +proof (subst integral_FTC_atLeastAtMost)
    1.21 +  fix x show "DERIV (\<lambda>x. x^Suc k / Suc k) x :> x^k"
    1.22 +    by (intro derivative_eq_intros) auto
    1.23 +qed (auto simp: field_simps)
    1.24 +
    1.25 +lemma has_bochner_integral_nn_integral:
    1.26 +  assumes "f \<in> borel_measurable M" "AE x in M. 0 \<le> f x"
    1.27 +  assumes "(\<integral>\<^sup>+x. f x \<partial>M) = ereal x"
    1.28 +  shows "has_bochner_integral M f x"
    1.29 +  unfolding has_bochner_integral_iff
    1.30 +proof
    1.31 +  show "integrable M f"
    1.32 +    using assms by (rule integrableI_nn_integral_finite)
    1.33 +qed (auto simp: assms integral_eq_nn_integral)
    1.34 +
    1.35 +lemma (in prob_space) distributed_AE2:
    1.36 +  assumes [measurable]: "distributed M N X f" "Measurable.pred N P"
    1.37 +  shows "(AE x in M. P (X x)) \<longleftrightarrow> (AE x in N. 0 < f x \<longrightarrow> P x)"
    1.38 +proof -
    1.39 +  have "(AE x in M. P (X x)) \<longleftrightarrow> (AE x in distr M N X. P x)"
    1.40 +    by (simp add: AE_distr_iff)
    1.41 +  also have "\<dots> \<longleftrightarrow> (AE x in density N f. P x)"
    1.42 +    unfolding distributed_distr_eq_density[OF assms(1)] ..
    1.43 +  also have "\<dots> \<longleftrightarrow>  (AE x in N. 0 < f x \<longrightarrow> P x)"
    1.44 +    by (rule AE_density) simp
    1.45 +  finally show ?thesis .
    1.46 +qed
    1.47 +
    1.48 +subsection {* Erlang *}
    1.49 +
    1.50 +lemma nn_intergal_power_times_exp_Icc:
    1.51 +  assumes [arith]: "0 \<le> a"
    1.52 +  shows "(\<integral>\<^sup>+x. ereal (x^k * exp (-x)) * indicator {0 .. a} x \<partial>lborel) =
    1.53 +    (1 - (\<Sum>n\<le>k. (a^n * exp (-a)) / fact n)) * fact k" (is "?I = _")
    1.54 +proof -
    1.55 +  let ?f = "\<lambda>k x. x^k * exp (-x) / fact k"
    1.56 +  let ?F = "\<lambda>k x. - (\<Sum>n\<le>k. (x^n * exp (-x)) / fact n)"
    1.57 +
    1.58 +  have "?I * (inverse (fact k)) = 
    1.59 +      (\<integral>\<^sup>+x. ereal (x^k * exp (-x)) * indicator {0 .. a} x * (inverse (fact k)) \<partial>lborel)"
    1.60 +    by (intro nn_integral_multc[symmetric]) auto
    1.61 +  also have "\<dots> = (\<integral>\<^sup>+x. ereal (?f k x) * indicator {0 .. a} x \<partial>lborel)"
    1.62 +    by (intro nn_integral_cong) (simp add: field_simps)
    1.63 +  also have "\<dots> = ereal (?F k a) - (?F k 0)"
    1.64 +  proof (rule nn_integral_FTC_atLeastAtMost)
    1.65 +    fix x assume "x \<in> {0..a}"
    1.66 +    show "DERIV (?F k) x :> ?f k x"
    1.67 +    proof(induction k)
    1.68 +      case 0 show ?case by (auto intro!: derivative_eq_intros)
    1.69 +    next
    1.70 +      case (Suc k)
    1.71 +      have "DERIV (\<lambda>x. ?F k x - (x^Suc k * exp (-x)) / fact (Suc k)) x :>
    1.72 +        ?f k x - ((real (Suc k) - x) * x ^ k * exp (- x)) / real (fact (Suc k))"
    1.73 +        by (intro DERIV_diff Suc)
    1.74 +           (auto intro!: derivative_eq_intros simp del: fact_Suc power_Suc
    1.75 +                 simp add: field_simps power_Suc[symmetric] real_of_nat_def[symmetric])
    1.76 +      also have "(\<lambda>x. ?F k x - (x^Suc k * exp (-x)) / fact (Suc k)) = ?F (Suc k)"
    1.77 +        by simp
    1.78 +      also have "?f k x - ((real (Suc k) - x) * x ^ k * exp (- x)) / real (fact (Suc k)) = ?f (Suc k) x"
    1.79 +        by (auto simp: field_simps simp del: fact_Suc)
    1.80 +           (simp_all add: real_of_nat_Suc field_simps)
    1.81 +      finally show ?case .
    1.82 +    qed
    1.83 +  qed auto
    1.84 +  also have "\<dots> = ereal (1 - (\<Sum>n\<le>k. (a^n * exp (-a)) / fact n))"
    1.85 +    by (auto simp: power_0_left if_distrib[where f="\<lambda>x. x / a" for a] setsum_cases)
    1.86 +  finally show ?thesis
    1.87 +    by (cases "?I") (auto simp: field_simps)
    1.88 +qed
    1.89 +
    1.90 +lemma nn_intergal_power_times_exp_Ici:
    1.91 +  shows "(\<integral>\<^sup>+x. ereal (x^k * exp (-x)) * indicator {0 ..} x \<partial>lborel) = fact k"
    1.92 +proof (rule LIMSEQ_unique)
    1.93 +  let ?X = "\<lambda>n. \<integral>\<^sup>+ x. ereal (x^k * exp (-x)) * indicator {0 .. real n} x \<partial>lborel"
    1.94 +  show "?X ----> (\<integral>\<^sup>+x. ereal (x^k * exp (-x)) * indicator {0 ..} x \<partial>lborel)"
    1.95 +    apply (intro nn_integral_LIMSEQ)
    1.96 +    apply (auto simp: incseq_def le_fun_def eventually_sequentially
    1.97 +                split: split_indicator intro!: Lim_eventually)
    1.98 +    apply (metis natceiling_le_eq)
    1.99 +    done
   1.100 +
   1.101 +  have "((\<lambda>x. (1 - (\<Sum>n\<le>k. (x ^ n / exp x) / real (fact n))) * fact k) ---> (1 - (\<Sum>n\<le>k. 0 / real (fact n))) * fact k) at_top"
   1.102 +    by (intro tendsto_intros tendsto_power_div_exp_0) simp
   1.103 +  then show "?X ----> fact k"
   1.104 +    by (subst nn_intergal_power_times_exp_Icc)
   1.105 +       (auto simp: exp_minus field_simps intro!: filterlim_compose[OF _ filterlim_real_sequentially])
   1.106 +qed
   1.107 +
   1.108 +definition erlang_density :: "nat \<Rightarrow> real \<Rightarrow> real \<Rightarrow> real" where
   1.109 +  "erlang_density k l x = (if x < 0 then 0 else (l^(Suc k) * x^k * exp (- l * x)) / fact k)"
   1.110 +
   1.111 +definition erlang_CDF ::  "nat \<Rightarrow> real \<Rightarrow> real \<Rightarrow> real" where
   1.112 +  "erlang_CDF k l x = (if x < 0 then 0 else 1 - (\<Sum>n\<le>k. ((l * x)^n * exp (- l * x) / fact n)))"
   1.113 +
   1.114 +lemma erlang_density_nonneg: "0 \<le> l \<Longrightarrow> 0 \<le> erlang_density k l x"
   1.115 +  by (simp add: erlang_density_def)
   1.116 +
   1.117 +lemma borel_measurable_erlang_density[measurable]: "erlang_density k l \<in> borel_measurable borel"
   1.118 +  by (auto simp add: erlang_density_def[abs_def])
   1.119 +
   1.120 +lemma erlang_CDF_transform: "0 < l \<Longrightarrow> erlang_CDF k l a = erlang_CDF k 1 (l * a)"
   1.121 +  by (auto simp add: erlang_CDF_def mult_less_0_iff)
   1.122 +
   1.123 +lemma nn_integral_erlang_density:
   1.124 +  assumes [arith]: "0 < l"
   1.125 +  shows "(\<integral>\<^sup>+ x. ereal (erlang_density k l x) * indicator {.. a} x \<partial>lborel) = erlang_CDF k l a"
   1.126 +proof cases
   1.127 +  assume [arith]: "0 \<le> a"
   1.128 +  have eq: "\<And>x. indicator {0..a} (x / l) = indicator {0..a*l} x"
   1.129 +    by (simp add: field_simps split: split_indicator)
   1.130 +  have "(\<integral>\<^sup>+x. ereal (erlang_density k l x) * indicator {.. a} x \<partial>lborel) =
   1.131 +    (\<integral>\<^sup>+x. (l/fact k) * (ereal ((l*x)^k * exp (- (l*x))) * indicator {0 .. a} x) \<partial>lborel)"
   1.132 +    by (intro nn_integral_cong) (auto simp: erlang_density_def power_mult_distrib split: split_indicator)
   1.133 +  also have "\<dots> = (l/fact k) * (\<integral>\<^sup>+x. ereal ((l*x)^k * exp (- (l*x))) * indicator {0 .. a} x \<partial>lborel)"
   1.134 +    by (intro nn_integral_cmult) auto
   1.135 +  also have "\<dots> = ereal (l/fact k) * ((1/l) * (\<integral>\<^sup>+x. ereal (x^k * exp (- x)) * indicator {0 .. l * a} x \<partial>lborel))"
   1.136 +    by (subst nn_integral_real_affine[where c="1 / l" and t=0]) (auto simp: field_simps eq)
   1.137 +  also have "\<dots> = (1 - (\<Sum>n\<le>k. ((l * a)^n * exp (-(l * a))) / fact n))"
   1.138 +    by (subst nn_intergal_power_times_exp_Icc) auto
   1.139 +  also have "\<dots> = erlang_CDF k l a"
   1.140 +    by (auto simp: erlang_CDF_def)
   1.141 +  finally show ?thesis .
   1.142 +next
   1.143 +  assume "\<not> 0 \<le> a" 
   1.144 +  moreover then have "(\<integral>\<^sup>+ x. ereal (erlang_density k l x) * indicator {.. a} x \<partial>lborel) = (\<integral>\<^sup>+x. 0 \<partial>(lborel::real measure))"
   1.145 +    by (intro nn_integral_cong) (auto simp: erlang_density_def)
   1.146 +  ultimately show ?thesis
   1.147 +    by (simp add: erlang_CDF_def)
   1.148 +qed
   1.149 +
   1.150 +lemma emeasure_erlang_density: 
   1.151 +  "0 < l \<Longrightarrow> emeasure (density lborel (erlang_density k l)) {.. a} = erlang_CDF k l a"
   1.152 +  by (simp add: emeasure_density nn_integral_erlang_density)
   1.153 +
   1.154 +lemma nn_integral_erlang_ith_moment: 
   1.155 +  fixes k i :: nat and l :: real
   1.156 +  assumes [arith]: "0 < l" 
   1.157 +  shows "(\<integral>\<^sup>+ x. ereal (erlang_density k l x * x ^ i) \<partial>lborel) = fact (k + i) / (fact k * l ^ i)"
   1.158 +proof -
   1.159 +  have eq: "\<And>x. indicator {0..} (x / l) = indicator {0..} x"
   1.160 +    by (simp add: field_simps split: split_indicator)
   1.161 +  have "(\<integral>\<^sup>+ x. ereal (erlang_density k l x * x^i) \<partial>lborel) =
   1.162 +    (\<integral>\<^sup>+x. (l/(fact k * l^i)) * (ereal ((l*x)^(k+i) * exp (- (l*x))) * indicator {0 ..} x) \<partial>lborel)"
   1.163 +    by (intro nn_integral_cong) (auto simp: erlang_density_def power_mult_distrib power_add split: split_indicator)
   1.164 +  also have "\<dots> = (l/(fact k * l^i)) * (\<integral>\<^sup>+x. ereal ((l*x)^(k+i) * exp (- (l*x))) * indicator {0 ..} x \<partial>lborel)"
   1.165 +    by (intro nn_integral_cmult) auto
   1.166 +  also have "\<dots> = ereal (l/(fact k * l^i)) * ((1/l) * (\<integral>\<^sup>+x. ereal (x^(k+i) * exp (- x)) * indicator {0 ..} x \<partial>lborel))"
   1.167 +    by (subst nn_integral_real_affine[where c="1 / l" and t=0]) (auto simp: field_simps eq)
   1.168 +  also have "\<dots> = fact (k + i) / (fact k * l ^ i)"
   1.169 +    by (subst nn_intergal_power_times_exp_Ici) auto
   1.170 +  finally show ?thesis .
   1.171 +qed
   1.172 +
   1.173 +lemma prob_space_erlang_density:
   1.174 +  assumes l[arith]: "0 < l"
   1.175 +  shows "prob_space (density lborel (erlang_density k l))" (is "prob_space ?D")
   1.176 +proof
   1.177 +  show "emeasure ?D (space ?D) = 1"
   1.178 +    using nn_integral_erlang_ith_moment[OF l, where k=k and i=0] by (simp add: emeasure_density)
   1.179 +qed
   1.180 +
   1.181 +lemma (in prob_space) erlang_distributed_le:
   1.182 +  assumes D: "distributed M lborel X (erlang_density k l)"
   1.183 +  assumes [simp, arith]: "0 < l" "0 \<le> a"
   1.184 +  shows "\<P>(x in M. X x \<le> a) = erlang_CDF k l a"
   1.185 +proof -
   1.186 +  have "emeasure M {x \<in> space M. X x \<le> a } = emeasure (distr M lborel X) {.. a}"
   1.187 +    using distributed_measurable[OF D]
   1.188 +    by (subst emeasure_distr) (auto intro!: arg_cong2[where f=emeasure])
   1.189 +  also have "\<dots> = emeasure (density lborel (erlang_density k l)) {.. a}"
   1.190 +    unfolding distributed_distr_eq_density[OF D] ..
   1.191 +  also have "\<dots> = erlang_CDF k l a"
   1.192 +    by (auto intro!: emeasure_erlang_density)
   1.193 +  finally show ?thesis
   1.194 +    by (auto simp: measure_def)
   1.195 +qed
   1.196 +
   1.197 +lemma (in prob_space) erlang_distributed_gt:
   1.198 +  assumes D[simp]: "distributed M lborel X (erlang_density k l)"
   1.199 +  assumes [arith]: "0 < l" "0 \<le> a"
   1.200 +  shows "\<P>(x in M. a < X x ) = 1 - (erlang_CDF k l a)"
   1.201 +proof -
   1.202 +  have " 1 - (erlang_CDF k l a) = 1 - \<P>(x in M. X x \<le> a)" by (subst erlang_distributed_le) auto
   1.203 +  also have "\<dots> = prob (space M - {x \<in> space M. X x \<le> a })"
   1.204 +    using distributed_measurable[OF D] by (auto simp: prob_compl)
   1.205 +  also have "\<dots> = \<P>(x in M. a < X x )" by (auto intro!: arg_cong[where f=prob] simp: not_le)
   1.206 +  finally show ?thesis by simp
   1.207 +qed
   1.208 +
   1.209 +lemma erlang_CDF_at0: "erlang_CDF k l 0 = 0"
   1.210 +  by (induction k) (auto simp: erlang_CDF_def)
   1.211 +
   1.212 +lemma erlang_distributedI:
   1.213 +  assumes X[measurable]: "X \<in> borel_measurable M" and [arith]: "0 < l"
   1.214 +    and X_distr: "\<And>a. 0 \<le> a \<Longrightarrow> emeasure M {x\<in>space M. X x \<le> a} = erlang_CDF k l a"
   1.215 +  shows "distributed M lborel X (erlang_density k l)"
   1.216 +proof (rule distributedI_borel_atMost)
   1.217 +  fix a :: real
   1.218 +  { assume "a \<le> 0"  
   1.219 +    with X have "emeasure M {x\<in>space M. X x \<le> a} \<le> emeasure M {x\<in>space M. X x \<le> 0}"
   1.220 +      by (intro emeasure_mono) auto
   1.221 +    also have "... = 0"  by (auto intro!: erlang_CDF_at0 simp: X_distr[of 0])
   1.222 +    finally have "emeasure M {x\<in>space M. X x \<le> a} \<le> 0" by simp
   1.223 +    then have "emeasure M {x\<in>space M. X x \<le> a} = 0" by (simp add:emeasure_le_0_iff)
   1.224 +  }
   1.225 +  note eq_0 = this
   1.226 +
   1.227 +  show "(\<integral>\<^sup>+ x. erlang_density k l x * indicator {..a} x \<partial>lborel) = ereal (erlang_CDF k l a)"
   1.228 +    using nn_integral_erlang_density[of l k a]
   1.229 +    by (simp add: times_ereal.simps(1)[symmetric] ereal_indicator del: times_ereal.simps)
   1.230 +
   1.231 +  show "emeasure M {x\<in>space M. X x \<le> a} = ereal (erlang_CDF k l a)"
   1.232 +    using X_distr[of a] eq_0 by (auto simp: one_ereal_def erlang_CDF_def)
   1.233 +qed (simp_all add: erlang_density_nonneg)
   1.234 +
   1.235 +lemma (in prob_space) erlang_distributed_iff:
   1.236 +  assumes [arith]: "0<l"
   1.237 +  shows "distributed M lborel X (erlang_density k l) \<longleftrightarrow>
   1.238 +    (X \<in> borel_measurable M \<and> 0 < l \<and>  (\<forall>a\<ge>0. \<P>(x in M. X x \<le> a) = erlang_CDF k l a ))"
   1.239 +  using
   1.240 +    distributed_measurable[of M lborel X "erlang_density k l"]
   1.241 +    emeasure_erlang_density[of l]
   1.242 +    erlang_distributed_le[of X k l]
   1.243 +  by (auto intro!: erlang_distributedI simp: one_ereal_def emeasure_eq_measure) 
   1.244 +
   1.245 +lemma (in prob_space) erlang_distributed_mult_const:
   1.246 +  assumes erlX: "distributed M lborel X (erlang_density k l)"
   1.247 +  assumes a_pos[arith]: "0 < \<alpha>"  "0 < l"
   1.248 +  shows  "distributed M lborel (\<lambda>x. \<alpha> * X x) (erlang_density k (l / \<alpha>))"
   1.249 +proof (subst erlang_distributed_iff, safe)
   1.250 +  have [measurable]: "random_variable borel X"  and  [arith]: "0 < l " 
   1.251 +  and  [simp]: "\<And>a. 0 \<le> a \<Longrightarrow> prob {x \<in> space M. X x \<le> a} = erlang_CDF k l a"
   1.252 +    by(insert erlX, auto simp: erlang_distributed_iff)
   1.253 +
   1.254 +  show "random_variable borel (\<lambda>x. \<alpha> * X x)" "0 < l / \<alpha>"  "0 < l / \<alpha>" 
   1.255 +    by (auto simp:field_simps)
   1.256 +  
   1.257 +  fix a:: real assume [arith]: "0 \<le> a"
   1.258 +  obtain b:: real  where [simp, arith]: "b = a/ \<alpha>" by blast 
   1.259 +
   1.260 +  have [arith]: "0 \<le> b" by (auto simp: divide_nonneg_pos)
   1.261 + 
   1.262 +  have "prob {x \<in> space M. \<alpha> * X x \<le> a}  = prob {x \<in> space M.  X x \<le> b}"
   1.263 +    by (rule arg_cong[where f= prob]) (auto simp:field_simps)
   1.264 +  
   1.265 +  moreover have "prob {x \<in> space M. X x \<le> b} = erlang_CDF k l b" by auto
   1.266 +  moreover have "erlang_CDF k (l / \<alpha>) a = erlang_CDF k l b" unfolding erlang_CDF_def by auto
   1.267 +  ultimately show "prob {x \<in> space M. \<alpha> * X x \<le> a} = erlang_CDF k (l / \<alpha>) a" by fastforce  
   1.268 +qed
   1.269 +
   1.270 +lemma (in prob_space) has_bochner_integral_erlang_ith_moment:
   1.271 +  fixes k i :: nat and l :: real
   1.272 +  assumes [arith]: "0 < l" and D: "distributed M lborel X (erlang_density k l)"
   1.273 +  shows "has_bochner_integral M (\<lambda>x. X x ^ i) (fact (k + i) / (fact k * l ^ i))"
   1.274 +proof (rule has_bochner_integral_nn_integral)
   1.275 +  show "AE x in M. 0 \<le> X x ^ i"
   1.276 +    by (subst distributed_AE2[OF D]) (auto simp: erlang_density_def)
   1.277 +  show "(\<integral>\<^sup>+ x. ereal (X x ^ i) \<partial>M) = ereal (fact (k + i) / (fact k * l ^ i))"
   1.278 +    using nn_integral_erlang_ith_moment[of l k i]
   1.279 +    by (subst distributed_nn_integral[symmetric, OF D]) auto
   1.280 +qed (insert distributed_measurable[OF D], simp)
   1.281 +
   1.282 +lemma (in prob_space) erlang_ith_moment_integrable:
   1.283 +  "0 < l \<Longrightarrow> distributed M lborel X (erlang_density k l) \<Longrightarrow> integrable M (\<lambda>x. X x ^ i)"
   1.284 +  by rule (rule has_bochner_integral_erlang_ith_moment)
   1.285 +
   1.286 +lemma (in prob_space) erlang_ith_moment:
   1.287 +  "0 < l \<Longrightarrow> distributed M lborel X (erlang_density k l) \<Longrightarrow>
   1.288 +    expectation (\<lambda>x. X x ^ i) = fact (k + i) / (fact k * l ^ i)"
   1.289 +  by (rule has_bochner_integral_integral_eq) (rule has_bochner_integral_erlang_ith_moment)
   1.290 +
   1.291 +lemma (in prob_space) erlang_distributed_variance:
   1.292 +  assumes [arith]: "0 < l" and "distributed M lborel X (erlang_density k l)"
   1.293 +  shows "variance X = (k + 1) / l\<^sup>2"
   1.294 +proof (subst variance_eq)
   1.295 +  show "integrable M X" "integrable M (\<lambda>x. (X x)\<^sup>2)"
   1.296 +    using erlang_ith_moment_integrable[OF assms, of 1] erlang_ith_moment_integrable[OF assms, of 2]
   1.297 +    by auto
   1.298 +
   1.299 +  show "expectation (\<lambda>x. (X x)\<^sup>2) - (expectation X)\<^sup>2 = real (k + 1) / l\<^sup>2"
   1.300 +    using erlang_ith_moment[OF assms, of 1] erlang_ith_moment[OF assms, of 2]
   1.301 +    by simp (auto simp: power2_eq_square field_simps real_of_nat_Suc)
   1.302 +qed
   1.303 +
   1.304  subsection {* Exponential distribution *}
   1.305  
   1.306 -definition exponential_density :: "real \<Rightarrow> real \<Rightarrow> real" where
   1.307 -  "exponential_density l x = (if x < 0 then 0 else l * exp (- x * l))"
   1.308 +abbreviation exponential_density :: "real \<Rightarrow> real \<Rightarrow> real" where
   1.309 +  "exponential_density \<equiv> erlang_density 0"
   1.310  
   1.311 -lemma borel_measurable_exponential_density[measurable]: "exponential_density l \<in> borel_measurable borel"
   1.312 -  by (auto simp add: exponential_density_def[abs_def])
   1.313 +lemma exponential_density_def:
   1.314 +  "exponential_density l x = (if x < 0 then 0 else l * exp (- x * l))"
   1.315 +  by (simp add: fun_eq_iff erlang_density_def)
   1.316 +
   1.317 +lemma erlang_CDF_0: "erlang_CDF 0 l a = (if 0 \<le> a then 1 - exp (- l * a) else 0)"
   1.318 +  by (simp add: erlang_CDF_def)
   1.319  
   1.320  lemma (in prob_space) exponential_distributed_params:
   1.321    assumes D: "distributed M lborel X (exponential_density l)"
   1.322 @@ -42,74 +341,20 @@
   1.323      by (simp add: emeasure_density distributed_distr_eq_density[OF D])
   1.324  qed assumption
   1.325  
   1.326 -lemma
   1.327 -  assumes [arith]: "0 < l"
   1.328 -  shows emeasure_exponential_density_le0: "0 \<le> a \<Longrightarrow> emeasure (density lborel (exponential_density l)) {.. a} = 1 - exp (- a * l)"
   1.329 -    and prob_space_exponential_density: "prob_space (density lborel (exponential_density l))"
   1.330 -      (is "prob_space ?D")
   1.331 -proof -
   1.332 -  let ?f = "\<lambda>x. l * exp (- x * l)"
   1.333 -  let ?F = "\<lambda>x. - exp (- x * l)"
   1.334 -
   1.335 -  have deriv: "\<And>x. DERIV ?F x :> ?f x" "\<And>x. 0 \<le> l * exp (- x * l)"
   1.336 -    by (auto intro!: derivative_eq_intros simp: zero_le_mult_iff)
   1.337 -
   1.338 -  have "emeasure ?D (space ?D) = (\<integral>\<^sup>+ x. ereal (?f x) * indicator {0..} x \<partial>lborel)"
   1.339 -    by (auto simp: emeasure_density exponential_density_def
   1.340 -             intro!: nn_integral_cong split: split_indicator)
   1.341 -  also have "\<dots> = ereal (0 - ?F 0)"
   1.342 -  proof (rule nn_integral_FTC_atLeast)
   1.343 -    have "((\<lambda>x. exp (l * x)) ---> 0) at_bot"
   1.344 -      by (rule filterlim_compose[OF exp_at_bot filterlim_tendsto_pos_mult_at_bot[of _ l]])
   1.345 -         (simp_all add: tendsto_const filterlim_ident)
   1.346 -    then show "((\<lambda>x. - exp (- x * l)) ---> 0) at_top"
   1.347 -      unfolding filterlim_at_top_mirror
   1.348 -      by (simp add: tendsto_minus_cancel_left[symmetric] ac_simps)
   1.349 -  qed (insert deriv, auto)
   1.350 -  also have "\<dots> = 1" by (simp add: one_ereal_def)
   1.351 -  finally have "emeasure ?D (space ?D) = 1" .
   1.352 -  then show "prob_space ?D" by rule
   1.353 -
   1.354 -  assume "0 \<le> a"
   1.355 -  have "emeasure ?D {..a} = (\<integral>\<^sup>+x. ereal (?f x) * indicator {0..a} x \<partial>lborel)"
   1.356 -    by (auto simp add: emeasure_density intro!: nn_integral_cong split: split_indicator)
   1.357 -       (auto simp: exponential_density_def)
   1.358 -  also have "(\<integral>\<^sup>+x. ereal (?f x) * indicator {0..a} x \<partial>lborel) = ereal (?F a) - ereal (?F 0)"
   1.359 -    using `0 \<le> a` deriv by (intro nn_integral_FTC_atLeastAtMost) auto
   1.360 -  also have "\<dots> = 1 - exp (- a * l)"
   1.361 -    by simp
   1.362 -  finally show "emeasure ?D {.. a} = 1 - exp (- a * l)" .
   1.363 -qed
   1.364 -
   1.365 +lemma prob_space_exponential_density: "0 < l \<Longrightarrow> prob_space (density lborel (exponential_density l))"
   1.366 +  by (rule prob_space_erlang_density)
   1.367  
   1.368  lemma (in prob_space) exponential_distributedD_le:
   1.369    assumes D: "distributed M lborel X (exponential_density l)" and a: "0 \<le> a"
   1.370    shows "\<P>(x in M. X x \<le> a) = 1 - exp (- a * l)"
   1.371 -proof -
   1.372 -  have "emeasure M {x \<in> space M. X x \<le> a } = emeasure (distr M lborel X) {.. a}"
   1.373 -    using distributed_measurable[OF D] 
   1.374 -    by (subst emeasure_distr) (auto intro!: arg_cong2[where f=emeasure])
   1.375 -  also have "\<dots> = emeasure (density lborel (exponential_density l)) {.. a}"
   1.376 -    unfolding distributed_distr_eq_density[OF D] ..
   1.377 -  also have "\<dots> = 1 - exp (- a * l)"
   1.378 -    using emeasure_exponential_density_le0[OF exponential_distributed_params[OF D] a] .
   1.379 -  finally show ?thesis
   1.380 -    by (auto simp: measure_def)
   1.381 -qed
   1.382 +  using erlang_distributed_le[OF D exponential_distributed_params[OF D] a] a
   1.383 +  by (simp add: erlang_CDF_def)
   1.384  
   1.385  lemma (in prob_space) exponential_distributedD_gt:
   1.386    assumes D: "distributed M lborel X (exponential_density l)" and a: "0 \<le> a"
   1.387    shows "\<P>(x in M. a < X x ) = exp (- a * l)"
   1.388 -proof -
   1.389 -  have "exp (- a * l) = 1 - \<P>(x in M. X x \<le> a)"
   1.390 -    unfolding exponential_distributedD_le[OF D a] by simp
   1.391 -  also have "\<dots> = prob (space M - {x \<in> space M. X x \<le> a })"
   1.392 -    using distributed_measurable[OF D]
   1.393 -    by (subst prob_compl) auto
   1.394 -  also have "\<dots> = \<P>(x in M. a < X x )"
   1.395 -    by (auto intro!: arg_cong[where f=prob] simp: not_le)
   1.396 -  finally show ?thesis by simp
   1.397 -qed
   1.398 +  using erlang_distributed_gt[OF D exponential_distributed_params[OF D] a] a
   1.399 +  by (simp add: erlang_CDF_def)
   1.400  
   1.401  lemma (in prob_space) exponential_distributed_memoryless:
   1.402    assumes D: "distributed M lborel X (exponential_density l)" and a: "0 \<le> a"and t: "0 \<le> t"
   1.403 @@ -129,100 +374,202 @@
   1.404    assumes X[measurable]: "X \<in> borel_measurable M" and [arith]: "0 < l"
   1.405      and X_distr: "\<And>a. 0 \<le> a \<Longrightarrow> emeasure M {x\<in>space M. X x \<le> a} = 1 - exp (- a * l)"
   1.406    shows "distributed M lborel X (exponential_density l)"
   1.407 -proof (rule distributedI_borel_atMost)
   1.408 -  fix a :: real
   1.409 -
   1.410 -  { assume "a \<le> 0"  
   1.411 -    with X have "emeasure M {x\<in>space M. X x \<le> a} \<le> emeasure M {x\<in>space M. X x \<le> 0}"
   1.412 -      by (intro emeasure_mono) auto
   1.413 -    then have "emeasure M {x\<in>space M. X x \<le> a} = 0"
   1.414 -      using X_distr[of 0] by (simp add: one_ereal_def emeasure_le_0_iff) }
   1.415 -  note eq_0 = this
   1.416 -
   1.417 -  have "\<And>x. \<not> 0 \<le> a \<Longrightarrow> ereal (exponential_density l x) * indicator {..a} x = 0"
   1.418 -    by (simp add: exponential_density_def)
   1.419 -  then show "(\<integral>\<^sup>+ x. exponential_density l x * indicator {..a} x \<partial>lborel) = ereal (if 0 \<le> a then 1 - exp (- a * l) else 0)"
   1.420 -    using emeasure_exponential_density_le0[of l a]
   1.421 -    by (auto simp: emeasure_density times_ereal.simps[symmetric] ereal_indicator
   1.422 -             simp del: times_ereal.simps ereal_zero_times)
   1.423 -  show "emeasure M {x\<in>space M. X x \<le> a} = ereal (if 0 \<le> a then 1 - exp (- a * l) else 0)"
   1.424 -    using X_distr[of a] eq_0 by (auto simp: one_ereal_def)
   1.425 -  show "AE x in lborel. 0 \<le> exponential_density l x "
   1.426 -    by (auto simp: exponential_density_def)
   1.427 -qed simp_all
   1.428 +proof (rule erlang_distributedI)
   1.429 +  fix a :: real assume "0 \<le> a" then show "emeasure M {x \<in> space M. X x \<le> a} = ereal (erlang_CDF 0 l a)"
   1.430 +    using X_distr[of a] by (simp add: erlang_CDF_def one_ereal_def)
   1.431 +qed fact+
   1.432  
   1.433  lemma (in prob_space) exponential_distributed_iff:
   1.434    "distributed M lborel X (exponential_density l) \<longleftrightarrow>
   1.435      (X \<in> borel_measurable M \<and> 0 < l \<and> (\<forall>a\<ge>0. \<P>(x in M. X x \<le> a) = 1 - exp (- a * l)))"
   1.436 -  using
   1.437 -    distributed_measurable[of M lborel X "exponential_density l"]
   1.438 -    exponential_distributed_params[of X l]
   1.439 -    emeasure_exponential_density_le0[of l]
   1.440 -    exponential_distributedD_le[of X l]
   1.441 -  by (auto intro!: exponential_distributedI simp: one_ereal_def emeasure_eq_measure)
   1.442 +  using exponential_distributed_params[of X l] erlang_distributed_iff[of l X 0] by (auto simp: erlang_CDF_0)
   1.443  
   1.444 -lemma borel_integral_x_exp:
   1.445 -  "has_bochner_integral lborel (\<lambda>x. x * exp (- x) * indicator {0::real ..} x) 1"
   1.446 -proof (rule has_bochner_integral_monotone_convergence)
   1.447 -  let ?f = "\<lambda>i x. x * exp (- x) * indicator {0::real .. i} x"
   1.448 -  have "eventually (\<lambda>b::real. 0 \<le> b) at_top"
   1.449 -    by (rule eventually_ge_at_top)
   1.450 -  then have "eventually (\<lambda>b. 1 - (inverse (exp b) + b / exp b) = integral\<^sup>L lborel (?f b)) at_top"
   1.451 -  proof eventually_elim
   1.452 -   fix b :: real assume [simp]: "0 \<le> b"
   1.453 -    have "(\<integral>x. (exp (-x)) * indicator {0 .. b} x \<partial>lborel) - (integral\<^sup>L lborel (?f b)) = 
   1.454 -      (\<integral>x. (exp (-x) - x * exp (-x)) * indicator {0 .. b} x \<partial>lborel)"
   1.455 -      by (subst integral_diff[symmetric])
   1.456 -         (auto intro!: borel_integrable_atLeastAtMost integral_cong split: split_indicator)
   1.457 -    also have "\<dots> = b * exp (-b) - 0 * exp (- 0)"
   1.458 -    proof (rule integral_FTC_atLeastAtMost)
   1.459 -      fix x assume "0 \<le> x" "x \<le> b"
   1.460 -      show "DERIV (\<lambda>x. x * exp (-x)) x :> exp (-x) - x * exp (-x)"
   1.461 -        by (auto intro!: derivative_eq_intros)
   1.462 -      show "isCont (\<lambda>x. exp (- x) - x * exp (- x)) x "
   1.463 -        by (intro continuous_intros)
   1.464 -    qed fact
   1.465 -    also have "(\<integral>x. (exp (-x)) * indicator {0 .. b} x \<partial>lborel) = - exp (- b) - - exp (- 0)"
   1.466 -      by (rule integral_FTC_atLeastAtMost) (auto intro!: derivative_eq_intros)
   1.467 -    finally show "1 - (inverse (exp b) + b / exp b) = integral\<^sup>L lborel (?f b)"
   1.468 -      by (auto simp: field_simps exp_minus inverse_eq_divide)
   1.469 -  qed
   1.470 -  then have "((\<lambda>i. integral\<^sup>L lborel (?f i)) ---> 1 - (0 + 0)) at_top"
   1.471 -  proof (rule Lim_transform_eventually)
   1.472 -    show "((\<lambda>i. 1 - (inverse (exp i) + i / exp i)) ---> 1 - (0 + 0 :: real)) at_top"
   1.473 -      using tendsto_power_div_exp_0[of 1]
   1.474 -      by (intro tendsto_intros tendsto_inverse_0_at_top exp_at_top) simp
   1.475 -  qed
   1.476 -  then show "(\<lambda>i. integral\<^sup>L lborel (?f i)) ----> 1"
   1.477 -    by (intro filterlim_compose[OF _ filterlim_real_sequentially]) simp
   1.478 -  show "AE x in lborel. mono (\<lambda>n::nat. x * exp (- x) * indicator {0..real n} x)"
   1.479 -    by (auto simp: mono_def mult_le_0_iff zero_le_mult_iff split: split_indicator)
   1.480 -  show "\<And>i::nat. integrable lborel (\<lambda>x. x * exp (- x) * indicator {0..real i} x)"
   1.481 -    by (rule borel_integrable_atLeastAtMost) auto
   1.482 -  show "AE x in lborel. (\<lambda>i. x * exp (- x) * indicator {0..real i} x) ----> x * exp (- x) * indicator {0..} x"
   1.483 -    apply (intro AE_I2 Lim_eventually )
   1.484 -    apply (rule_tac c="natfloor x + 1" in eventually_sequentiallyI)
   1.485 -    apply (auto simp add: not_le dest!: ge_natfloor_plus_one_imp_gt[simplified] split: split_indicator)
   1.486 -    done
   1.487 -qed (auto intro!: borel_measurable_times borel_measurable_exp)
   1.488  
   1.489  lemma (in prob_space) exponential_distributed_expectation:
   1.490 -  assumes D: "distributed M lborel X (exponential_density l)"
   1.491 -  shows "expectation X = 1 / l"
   1.492 -proof (subst distributed_integral[OF D, of "\<lambda>x. x", symmetric])
   1.493 -  have "0 < l"
   1.494 -   using exponential_distributed_params[OF D] .
   1.495 -  have [simp]: "\<And>x. x * (l * (exp (- (x * l)) * indicator {0..} (x * l))) =
   1.496 -    x * exponential_density l x"
   1.497 -    using `0 < l`
   1.498 -    by (auto split: split_indicator simp: zero_le_mult_iff exponential_density_def)
   1.499 -  from borel_integral_x_exp `0 < l`
   1.500 -  have "has_bochner_integral lborel (\<lambda>x. exponential_density l x * x) (1 / l)"
   1.501 -    by (subst (asm) lborel_has_bochner_integral_real_affine_iff[of l _ _ 0])
   1.502 -       (simp_all add: field_simps)
   1.503 -  then show "(\<integral> x. exponential_density l x * x \<partial>lborel) = 1 / l"
   1.504 -    by (metis has_bochner_integral_integral_eq)
   1.505 -qed simp
   1.506 +  "distributed M lborel X (exponential_density l) \<Longrightarrow> expectation X = 1 / l"
   1.507 +  using erlang_ith_moment[OF exponential_distributed_params, of X l X 0 1] by simp
   1.508 +
   1.509 +lemma exponential_density_nonneg: "0 < l \<Longrightarrow> 0 \<le> exponential_density l x"
   1.510 +  by (auto simp: exponential_density_def)
   1.511 +
   1.512 +lemma (in prob_space) exponential_distributed_min:
   1.513 +  assumes expX: "distributed M lborel X (exponential_density l)"
   1.514 +  assumes expY: "distributed M lborel Y (exponential_density u)"
   1.515 +  assumes ind: "indep_var borel X borel Y"
   1.516 +  shows "distributed M lborel (\<lambda>x. min (X x) (Y x)) (exponential_density (l + u))"
   1.517 +proof (subst exponential_distributed_iff, safe)
   1.518 +  have randX: "random_variable borel X" using expX by (simp add: exponential_distributed_iff)
   1.519 +  moreover have randY: "random_variable borel Y" using expY by (simp add: exponential_distributed_iff)
   1.520 +  ultimately show "random_variable borel (\<lambda>x. min (X x) (Y x))" by auto
   1.521 +  
   1.522 +  have "0 < l" by (rule exponential_distributed_params) fact
   1.523 +  moreover have "0 < u" by (rule exponential_distributed_params) fact
   1.524 +  ultimately  show " 0 < l + u" by force
   1.525 +
   1.526 +  fix a::real assume a[arith]: "0 \<le> a"
   1.527 +  have gt1[simp]: "\<P>(x in M. a < X x ) = exp (- a * l)" by (rule exponential_distributedD_gt[OF expX a]) 
   1.528 +  have gt2[simp]: "\<P>(x in M. a < Y x ) = exp (- a * u)" by (rule exponential_distributedD_gt[OF expY a]) 
   1.529 +
   1.530 +  have "\<P>(x in M. a < (min (X x) (Y x)) ) =  \<P>(x in M. a < (X x) \<and> a < (Y x))" by (auto intro!:arg_cong[where f=prob])
   1.531 +
   1.532 +  also have " ... =  \<P>(x in M. a < (X x)) *  \<P>(x in M. a< (Y x) )"
   1.533 +    using prob_indep_random_variable[OF ind, of "{a <..}" "{a <..}"] by simp
   1.534 +  also have " ... = exp (- a * (l + u))" by (auto simp:field_simps mult_exp_exp)
   1.535 +  finally have indep_prob: "\<P>(x in M. a < (min (X x) (Y x)) ) = exp (- a * (l + u))" .
   1.536 +
   1.537 +  have "{x \<in> space M. (min (X x) (Y x)) \<le>a } = (space M - {x \<in> space M. a<(min (X x) (Y x)) })"
   1.538 +    by auto
   1.539 +  then have "1 - prob {x \<in> space M. a < (min (X x) (Y x))} = prob {x \<in> space M. (min (X x) (Y x)) \<le> a}"
   1.540 +    using randX randY by (auto simp: prob_compl) 
   1.541 +  then show "prob {x \<in> space M. (min (X x) (Y x)) \<le> a} = 1 - exp (- a * (l + u))"
   1.542 +    using indep_prob by auto
   1.543 +qed
   1.544 + 
   1.545 +lemma (in prob_space) exponential_distributed_Min:
   1.546 +  assumes finI: "finite I"
   1.547 +  assumes A: "I \<noteq> {}"
   1.548 +  assumes expX: "\<And>i. i \<in> I \<Longrightarrow> distributed M lborel (X i) (exponential_density (l i))"
   1.549 +  assumes ind: "indep_vars (\<lambda>i. borel) X I" 
   1.550 +  shows "distributed M lborel (\<lambda>x. Min ((\<lambda>i. X i x)`I)) (exponential_density (\<Sum>i\<in>I. l i))"
   1.551 +using assms
   1.552 +proof (induct rule: finite_ne_induct)
   1.553 +  case (singleton i) then show ?case by simp
   1.554 +next
   1.555 +  case (insert i I)
   1.556 +  then have "distributed M lborel (\<lambda>x. min (X i x) (Min ((\<lambda>i. X i x)`I))) (exponential_density (l i + (\<Sum>i\<in>I. l i)))"
   1.557 +      by (intro exponential_distributed_min indep_vars_Min insert)
   1.558 +         (auto intro: indep_vars_subset) 
   1.559 +  then show ?case
   1.560 +    using insert by simp
   1.561 +qed
   1.562 +
   1.563 +lemma (in prob_space) exponential_distributed_variance:
   1.564 +  "distributed M lborel X (exponential_density l) \<Longrightarrow> variance X = 1 / l\<^sup>2"
   1.565 +  using erlang_distributed_variance[OF exponential_distributed_params, of X l X 0] by simp
   1.566 +
   1.567 +lemma nn_integral_zero': "AE x in M. f x = 0 \<Longrightarrow> (\<integral>\<^sup>+x. f x \<partial>M) = 0"
   1.568 +  by (simp cong: nn_integral_cong_AE)
   1.569 +
   1.570 +lemma convolution_erlang_density:
   1.571 +  fixes k\<^sub>1 k\<^sub>2 :: nat
   1.572 +  assumes [simp, arith]: "0 < l"
   1.573 +  shows "(\<lambda>x. \<integral>\<^sup>+y. ereal (erlang_density k\<^sub>1 l (x - y)) * ereal (erlang_density k\<^sub>2 l y) \<partial>lborel) =
   1.574 +    (erlang_density (Suc k\<^sub>1 + Suc k\<^sub>2 - 1) l)"
   1.575 +      (is "?LHS = ?RHS")
   1.576 +proof
   1.577 +  fix x :: real
   1.578 +  have "x \<le> 0 \<or> 0 < x"
   1.579 +    by arith
   1.580 +  then show "?LHS x = ?RHS x"
   1.581 +  proof
   1.582 +    assume "x \<le> 0" then show ?thesis
   1.583 +      apply (subst nn_integral_zero')
   1.584 +      apply (rule AE_I[where N="{0}"])
   1.585 +      apply (auto simp add: erlang_density_def not_less)
   1.586 +      done
   1.587 +  next
   1.588 +    note zero_le_mult_iff[simp] zero_le_divide_iff[simp]
   1.589 +  
   1.590 +    have I_eq1: "integral\<^sup>N lborel (erlang_density (Suc k\<^sub>1 + Suc k\<^sub>2 - 1) l) = 1"
   1.591 +      using nn_integral_erlang_ith_moment[of l "Suc k\<^sub>1 + Suc k\<^sub>2 - 1" 0] by (simp del: fact_Suc)
   1.592 +  
   1.593 +    have 1: "(\<integral>\<^sup>+ x. ereal (erlang_density (Suc k\<^sub>1 + Suc k\<^sub>2 - 1) l x * indicator {0<..} x) \<partial>lborel) = 1"
   1.594 +      apply (subst I_eq1[symmetric])
   1.595 +      unfolding erlang_density_def
   1.596 +      by (auto intro!: nn_integral_cong split:split_indicator)
   1.597 +  
   1.598 +    have "prob_space (density lborel ?LHS)"
   1.599 +      unfolding times_ereal.simps[symmetric]
   1.600 +      by (intro prob_space_convolution_density) 
   1.601 +         (auto intro!: prob_space_erlang_density erlang_density_nonneg)
   1.602 +    then have 2: "integral\<^sup>N lborel ?LHS = 1"
   1.603 +      by (auto dest!: prob_space.emeasure_space_1 simp: emeasure_density)
   1.604 +  
   1.605 +    let ?I = "(integral\<^sup>N lborel (\<lambda>y. ereal ((1 - y)^ k\<^sub>1 * y^k\<^sub>2 * indicator {0..1} y)))"
   1.606 +    let ?C = "real (fact (Suc (k\<^sub>1 + k\<^sub>2))) / (real (fact k\<^sub>1) * real (fact k\<^sub>2))"
   1.607 +    let ?s = "Suc k\<^sub>1 + Suc k\<^sub>2 - 1"
   1.608 +    let ?L = "(\<lambda>x. \<integral>\<^sup>+y. ereal (erlang_density k\<^sub>1 l (x- y) * erlang_density k\<^sub>2 l y * indicator {0..x} y) \<partial>lborel)"
   1.609 +
   1.610 +    { fix x :: real assume [arith]: "0 < x"
   1.611 +      have *: "\<And>x y n. (x - y * x::real)^n = x^n * (1 - y)^n"
   1.612 +        unfolding power_mult_distrib[symmetric] by (simp add: field_simps)
   1.613 +    
   1.614 +      have "?LHS x = ?L x"
   1.615 +        unfolding erlang_density_def
   1.616 +        by (auto intro!: nn_integral_cong split:split_indicator)
   1.617 +      also have "... = (\<lambda>x. ereal ?C * ?I * erlang_density ?s l x) x"
   1.618 +        apply (subst nn_integral_real_affine[where c=x and t = 0])
   1.619 +        apply (simp_all add: nn_integral_cmult[symmetric] nn_integral_multc[symmetric] erlang_density_nonneg del: fact_Suc)
   1.620 +        apply (intro nn_integral_cong)
   1.621 +        apply (auto simp add: erlang_density_def mult_less_0_iff exp_minus field_simps exp_diff power_add *
   1.622 +                    simp del: fact_Suc split: split_indicator)
   1.623 +        done
   1.624 +      finally have "(\<integral>\<^sup>+y. ereal (erlang_density k\<^sub>1 l (x - y) * erlang_density k\<^sub>2 l y) \<partial>lborel) = 
   1.625 +        (\<lambda>x. ereal ?C * ?I * erlang_density ?s l x) x"
   1.626 +        by simp }
   1.627 +    note * = this
   1.628 +
   1.629 +    assume [arith]: "0 < x"
   1.630 +    have 3: "1 = integral\<^sup>N lborel (\<lambda>xa. ?LHS xa * indicator {0<..} xa)"
   1.631 +      by (subst 2[symmetric])
   1.632 +         (auto intro!: nn_integral_cong_AE AE_I[where N="{0}"]
   1.633 +               simp: erlang_density_def  nn_integral_multc[symmetric] indicator_def split: split_if_asm)
   1.634 +    also have "... = integral\<^sup>N lborel (\<lambda>x. (ereal (?C) * ?I) * ((erlang_density ?s l x) * indicator {0<..} x))"
   1.635 +      by (auto intro!: nn_integral_cong simp: * split: split_indicator)
   1.636 +    also have "... = ereal (?C) * ?I"
   1.637 +      using 1
   1.638 +      by (auto simp: nn_integral_nonneg nn_integral_cmult)  
   1.639 +    finally have " ereal (?C) * ?I = 1" by presburger
   1.640 +  
   1.641 +    then show ?thesis
   1.642 +      using * by simp
   1.643 +  qed
   1.644 +qed
   1.645 +
   1.646 +lemma (in prob_space) sum_indep_erlang:
   1.647 +  assumes indep: "indep_var borel X borel Y"
   1.648 +  assumes [simp, arith]: "0 < l"
   1.649 +  assumes erlX: "distributed M lborel X (erlang_density k\<^sub>1 l)"
   1.650 +  assumes erlY: "distributed M lborel Y (erlang_density k\<^sub>2 l)"
   1.651 +  shows "distributed M lborel (\<lambda>x. X x + Y x) (erlang_density (Suc k\<^sub>1 + Suc k\<^sub>2 - 1) l)"
   1.652 +  using assms
   1.653 +  apply (subst convolution_erlang_density[symmetric, OF `0<l`])
   1.654 +  apply (intro distributed_convolution)
   1.655 +  apply auto
   1.656 +  done
   1.657 +
   1.658 +lemma (in prob_space) erlang_distributed_setsum:
   1.659 +  assumes finI : "finite I"
   1.660 +  assumes A: "I \<noteq> {}"
   1.661 +  assumes [simp, arith]: "0 < l"
   1.662 +  assumes expX: "\<And>i. i \<in> I \<Longrightarrow> distributed M lborel (X i) (erlang_density (k i) l)"
   1.663 +  assumes ind: "indep_vars (\<lambda>i. borel) X I"
   1.664 +  shows "distributed M lborel (\<lambda>x. \<Sum>i\<in>I. X i x) (erlang_density ((\<Sum>i\<in>I. Suc (k i)) - 1) l)"
   1.665 +using assms
   1.666 +proof (induct rule: finite_ne_induct)
   1.667 +  case (singleton i) then show ?case by auto
   1.668 +next
   1.669 +  case (insert i I)
   1.670 +    then have "distributed M lborel (\<lambda>x. (X i x) + (\<Sum>i\<in> I. X i x)) (erlang_density (Suc (k i) + Suc ((\<Sum>i\<in>I. Suc (k i)) - 1) - 1) l)"
   1.671 +      by(intro sum_indep_erlang indep_vars_setsum) (auto intro!: indep_vars_subset)
   1.672 +    also have "(\<lambda>x. (X i x) + (\<Sum>i\<in> I. X i x)) = (\<lambda>x. \<Sum>i\<in>insert i I. X i x)"
   1.673 +      using insert by auto
   1.674 +    also have "Suc(k i) + Suc ((\<Sum>i\<in>I. Suc (k i)) - 1) - 1 = (\<Sum>i\<in>insert i I. Suc (k i)) - 1"
   1.675 +      using insert by (auto intro!: Suc_pred simp: ac_simps)    
   1.676 +    finally show ?case by fast
   1.677 +qed
   1.678 +
   1.679 +lemma (in prob_space) exponential_distributed_setsum:
   1.680 +  assumes finI: "finite I"
   1.681 +  assumes A: "I \<noteq> {}"
   1.682 +  assumes expX: "\<And>i. i \<in> I \<Longrightarrow> distributed M lborel (X i) (exponential_density l)"
   1.683 +  assumes ind: "indep_vars (\<lambda>i. borel) X I" 
   1.684 +  shows "distributed M lborel (\<lambda>x. \<Sum>i\<in>I. X i x) (erlang_density ((card I) - 1) l)"
   1.685 +proof -
   1.686 +  obtain i where "i \<in> I" using assms by auto
   1.687 +  note exponential_distributed_params[OF expX[OF this]]
   1.688 +  from erlang_distributed_setsum[OF assms(1,2) this assms(3,4)] show ?thesis by simp
   1.689 +qed
   1.690  
   1.691  subsection {* Uniform distribution *}
   1.692  
   1.693 @@ -393,4 +740,19 @@
   1.694    finally show "(\<integral> x. indicator {a .. b} x / measure lborel {a .. b} * x \<partial>lborel) = (a + b) / 2" .
   1.695  qed auto
   1.696  
   1.697 +lemma (in prob_space) uniform_distributed_variance:
   1.698 +  fixes a b :: real
   1.699 +  assumes D: "distributed M lborel X (\<lambda>x. indicator {a .. b} x / measure lborel {a .. b})"
   1.700 +  shows "variance X = (b - a)\<^sup>2 / 12"
   1.701 +proof (subst distributed_variance)
   1.702 +  have [arith]: "a < b" using uniform_distributed_bounds[OF D] .
   1.703 +  let ?\<mu> = "expectation X" let ?D = "\<lambda>x. indicator {a..b} (x + ?\<mu>) / measure lborel {a..b}"
   1.704 +  have "(\<integral>x. x\<^sup>2 * (?D x) \<partial>lborel) = (\<integral>x. x\<^sup>2 * (indicator {a - ?\<mu> .. b - ?\<mu>} x) / measure lborel {a .. b} \<partial>lborel)"
   1.705 +    by (intro integral_cong) (auto split: split_indicator)
   1.706 +  also have "\<dots> = (b - a)\<^sup>2 / 12"
   1.707 +    by (simp add: integral_power measure_lebesgue_Icc uniform_distributed_expectation[OF D])
   1.708 +       (simp add: eval_nat_numeral field_simps )
   1.709 +  finally show "(\<integral>x. x\<^sup>2 * ?D x \<partial>lborel) = (b - a)\<^sup>2 / 12" .
   1.710 +qed fact
   1.711 +
   1.712  end