lemmas about the moments of the normal distribution
authorhoelzl
Mon Jun 16 13:19:48 2014 +0200 (2014-06-16)
changeset 57254d3d91422f408
parent 57253 6515cf25de13
child 57255 488046fdda59
lemmas about the moments of the normal distribution
CONTRIBUTORS
src/HOL/Probability/Distributions.thy
     1.1 --- a/CONTRIBUTORS	Fri Jun 13 14:49:59 2014 +0100
     1.2 +++ b/CONTRIBUTORS	Mon Jun 16 13:19:48 2014 +0200
     1.3 @@ -11,7 +11,7 @@
     1.4    Work on exotic automatic theorem provers for Sledgehammer (LEO-II, veriT,
     1.5    Waldmeister, etc.).
     1.6  
     1.7 -* June 2014: Sudeep Kanav, TUM, and Johannes Hölzl, TUM
     1.8 +* June 2014: Sudeep Kanav, TUM, Jeremy Avigad, CMU, and Johannes Hölzl, TUM
     1.9    Various properties of exponentially, Erlang, and normal distributed random variables.
    1.10  
    1.11  * May 2014: Cezary Kaliszyk, University of Innsbruck, and Jasmin Blanchette, TUM
     2.1 --- a/src/HOL/Probability/Distributions.thy	Fri Jun 13 14:49:59 2014 +0100
     2.2 +++ b/src/HOL/Probability/Distributions.thy	Mon Jun 16 13:19:48 2014 +0200
     2.3 @@ -1,6 +1,7 @@
     2.4  (*  Title:      HOL/Probability/Distributions.thy
     2.5      Author:     Sudeep Kanav, TU München
     2.6 -    Author:     Johannes Hölzl, TU München *)
     2.7 +    Author:     Johannes Hölzl, TU München
     2.8 +    Author:     Jeremy Avigad, CMU *)
     2.9  
    2.10  header {* Properties of Various Distributions *}
    2.11  
    2.12 @@ -8,36 +9,224 @@
    2.13    imports Convolution Information
    2.14  begin
    2.15  
    2.16 -lemma nn_integral_even_function:
    2.17 -  fixes f :: "real \<Rightarrow> ereal"
    2.18 -  assumes [measurable]: "f \<in> borel_measurable borel"
    2.19 -  assumes even: "\<And>x. f x = f (- x)"
    2.20 -  shows "(\<integral>\<^sup>+x. f x \<partial>lborel) = 2 * (\<integral>\<^sup>+x. f x * indicator {0 ..} x \<partial>lborel)"
    2.21 -proof -
    2.22 -  def f' \<equiv> "\<lambda>x. max 0 (f x)"
    2.23 -  have [measurable]: "f' \<in> borel_measurable borel"
    2.24 -    by (simp add: f'_def)
    2.25 +lemma tendsto_at_topI_sequentially:
    2.26 +  fixes f :: "real \<Rightarrow> 'b::first_countable_topology"
    2.27 +  assumes *: "\<And>X. filterlim X at_top sequentially \<Longrightarrow> (\<lambda>n. f (X n)) ----> y"
    2.28 +  shows "(f ---> y) at_top"
    2.29 +  unfolding filterlim_iff
    2.30 +proof safe
    2.31 +  fix P assume "eventually P (nhds y)"
    2.32 +  then have seq: "\<And>f. f ----> y \<Longrightarrow> eventually (\<lambda>x. P (f x)) at_top"
    2.33 +    unfolding eventually_nhds_iff_sequentially by simp
    2.34 +
    2.35 +  show "eventually (\<lambda>x. P (f x)) at_top"
    2.36 +  proof (rule ccontr)
    2.37 +    assume "\<not> eventually (\<lambda>x. P (f x)) at_top"
    2.38 +    then have "\<And>X. \<exists>x>X. \<not> P (f x)"
    2.39 +      unfolding eventually_at_top_dense by simp
    2.40 +    then obtain r where not_P: "\<And>x. \<not> P (f (r x))" and r: "\<And>x. x < r x"
    2.41 +      by metis
    2.42 +    
    2.43 +    def s \<equiv> "rec_nat (r 0) (\<lambda>_ x. r (x + 1))"
    2.44 +    then have [simp]: "s 0 = r 0" "\<And>n. s (Suc n) = r (s n + 1)"
    2.45 +      by auto
    2.46 +
    2.47 +    { fix n have "n < s n" using r
    2.48 +        by (induct n) (auto simp add: real_of_nat_Suc intro: less_trans add_strict_right_mono) }
    2.49 +    note s_subseq = this
    2.50 +
    2.51 +    have "mono s"
    2.52 +    proof (rule incseq_SucI)
    2.53 +      fix n show "s n \<le> s (Suc n)"
    2.54 +        using r[of "s n + 1"] by simp
    2.55 +    qed
    2.56 +
    2.57 +    have "filterlim s at_top sequentially"
    2.58 +      unfolding filterlim_at_top_gt[where c=0] eventually_sequentially
    2.59 +    proof (safe intro!: exI)
    2.60 +      fix Z :: real and n assume "0 < Z" "natceiling Z \<le> n"
    2.61 +      with real_natceiling_ge[of Z] `n < s n`
    2.62 +      show "Z \<le> s n"
    2.63 +        by auto
    2.64 +    qed
    2.65 +    moreover then have "eventually (\<lambda>x. P (f (s x))) sequentially"
    2.66 +      by (rule seq[OF *])
    2.67 +    then obtain n where "P (f (s n))"
    2.68 +      by (auto simp: eventually_sequentially)
    2.69 +    then show False
    2.70 +      using not_P by (cases n) auto
    2.71 +  qed
    2.72 +qed
    2.73    
    2.74 -  { fix x :: ereal have "2 * x = x + x"
    2.75 -      by (cases x) auto }
    2.76 -  note two_mult = this
    2.77 +lemma tendsto_integral_at_top:
    2.78 +  fixes f :: "real \<Rightarrow> 'a::{banach, second_countable_topology}"
    2.79 +  assumes [simp]: "sets M = sets borel" and f[measurable]: "integrable M f"
    2.80 +  shows "((\<lambda>y. \<integral> x. indicator {.. y} x *\<^sub>R f x \<partial>M) ---> \<integral> x. f x \<partial>M) at_top"
    2.81 +proof (rule tendsto_at_topI_sequentially)
    2.82 +  fix X :: "nat \<Rightarrow> real" assume "filterlim X at_top sequentially"
    2.83 +  show "(\<lambda>n. \<integral>x. indicator {..X n} x *\<^sub>R f x \<partial>M) ----> integral\<^sup>L M f"
    2.84 +  proof (rule integral_dominated_convergence)
    2.85 +    show "integrable M (\<lambda>x. norm (f x))"
    2.86 +      by (rule integrable_norm) fact
    2.87 +    show "AE x in M. (\<lambda>n. indicator {..X n} x *\<^sub>R f x) ----> f x"
    2.88 +    proof
    2.89 +      fix x
    2.90 +      from `filterlim X at_top sequentially` 
    2.91 +      have "eventually (\<lambda>n. x \<le> X n) sequentially"
    2.92 +        unfolding filterlim_at_top_ge[where c=x] by auto
    2.93 +      then show "(\<lambda>n. indicator {..X n} x *\<^sub>R f x) ----> f x"
    2.94 +        by (intro Lim_eventually) (auto split: split_indicator elim!: eventually_elim1)
    2.95 +    qed
    2.96 +    fix n show "AE x in M. norm (indicator {..X n} x *\<^sub>R f x) \<le> norm (f x)"
    2.97 +      by (auto split: split_indicator)
    2.98 +  qed auto
    2.99 +qed
   2.100 +
   2.101 +lemma filterlim_at_top_imp_at_infinity:
   2.102 +  fixes f :: "_ \<Rightarrow> real"
   2.103 +  shows "filterlim f at_top F \<Longrightarrow> filterlim f at_infinity F"
   2.104 +  by (rule filterlim_mono[OF _ at_top_le_at_infinity order_refl])
   2.105 +
   2.106 +lemma measurable_discrete_difference:
   2.107 +  fixes f :: "'a \<Rightarrow> 'b::t1_space"
   2.108 +  assumes f: "f \<in> measurable M N"
   2.109 +  assumes X: "countable X"
   2.110 +  assumes sets: "\<And>x. x \<in> X \<Longrightarrow> {x} \<in> sets M"
   2.111 +  assumes space: "\<And>x. x \<in> X \<Longrightarrow> g x \<in> space N"
   2.112 +  assumes eq: "\<And>x. x \<in> space M \<Longrightarrow> x \<notin> X \<Longrightarrow> f x = g x"
   2.113 +  shows "g \<in> measurable M N"
   2.114 +  unfolding measurable_def
   2.115 +proof safe
   2.116 +  fix x assume "x \<in> space M" then show "g x \<in> space N"
   2.117 +    using measurable_space[OF f, of x] eq[of x] space[of x] by (cases "x \<in> X") auto
   2.118 +next
   2.119 +  fix S assume S: "S \<in> sets N"
   2.120 +  have "g -` S \<inter> space M = (f -` S \<inter> space M) - (\<Union>x\<in>X. {x}) \<union> (\<Union>x\<in>{x\<in>X. g x \<in> S}. {x})"
   2.121 +    using sets.sets_into_space[OF sets] eq by auto
   2.122 +  also have "\<dots> \<in> sets M"
   2.123 +    by (safe intro!: sets.Diff sets.Un measurable_sets[OF f] S sets.countable_UN' X countable_Collect sets)
   2.124 +  finally show "g -` S \<inter> space M \<in> sets M" .
   2.125 +qed
   2.126  
   2.127 -  have "(\<integral>\<^sup>+x. f x \<partial>lborel) = (\<integral>\<^sup>+x. f' x \<partial>lborel)"
   2.128 -    unfolding f'_def nn_integral_max_0 ..
   2.129 -  also have "\<dots> = (\<integral>\<^sup>+x. f' x * indicator {0 ..} x + f' x * indicator {.. 0} x \<partial>lborel)"
   2.130 -    by (intro nn_integral_cong_AE AE_I[where N="{0}"]) (auto split: split_indicator_asm)
   2.131 -  also have "\<dots> = (\<integral>\<^sup>+x. f' x * indicator {0 ..} x \<partial>lborel) + (\<integral>\<^sup>+x. f' x * indicator {.. 0} x \<partial>lborel)"
   2.132 -    by (intro nn_integral_add) (auto simp: f'_def)
   2.133 -  also have "(\<integral>\<^sup>+x. f' x * indicator {.. 0} x \<partial>lborel) = (\<integral>\<^sup>+x. f' x * indicator {0 ..} x \<partial>lborel)"
   2.134 -    using even
   2.135 -    by (subst nn_integral_real_affine[where c="-1" and t=0])
   2.136 -       (auto simp: f'_def one_ereal_def[symmetric] split: split_indicator intro!: nn_integral_cong)
   2.137 -  also have "(\<integral>\<^sup>+x. f' x * indicator {0 ..} x \<partial>lborel) = (\<integral>\<^sup>+x. f x * indicator {0 ..} x \<partial>lborel)"
   2.138 -    unfolding f'_def by (subst (2) nn_integral_max_0[symmetric]) (auto intro!: nn_integral_cong split: split_indicator split_max)
   2.139 -  finally show ?thesis
   2.140 -    unfolding two_mult by simp
   2.141 +lemma AE_discrete_difference:
   2.142 +  assumes X: "countable X"
   2.143 +  assumes null: "\<And>x. x \<in> X \<Longrightarrow> emeasure M {x} = 0" 
   2.144 +  assumes sets: "\<And>x. x \<in> X \<Longrightarrow> {x} \<in> sets M"
   2.145 +  shows "AE x in M. x \<notin> X"
   2.146 +proof -
   2.147 +  have X_sets: "(\<Union>x\<in>X. {x}) \<in> sets M"
   2.148 +    using assms by (intro sets.countable_UN') auto
   2.149 +  have "emeasure M (\<Union>x\<in>X. {x}) = (\<integral>\<^sup>+ i. emeasure M {i} \<partial>count_space X)"
   2.150 +    by (rule emeasure_UN_countable) (auto simp: assms disjoint_family_on_def)
   2.151 +  also have "\<dots> = (\<integral>\<^sup>+ i. 0 \<partial>count_space X)"
   2.152 +    by (intro nn_integral_cong) (simp add: null)
   2.153 +  finally show "AE x in M. x \<notin> X"
   2.154 +    using AE_iff_measurable[of X M "\<lambda>x. x \<notin> X"] X_sets sets.sets_into_space[OF sets] by auto
   2.155 +qed
   2.156 +
   2.157 +lemma integrable_discrete_difference:
   2.158 +  fixes f :: "'a \<Rightarrow> 'b::{banach, second_countable_topology}"
   2.159 +  assumes X: "countable X"
   2.160 +  assumes null: "\<And>x. x \<in> X \<Longrightarrow> emeasure M {x} = 0" 
   2.161 +  assumes sets: "\<And>x. x \<in> X \<Longrightarrow> {x} \<in> sets M"
   2.162 +  assumes eq: "\<And>x. x \<in> space M \<Longrightarrow> x \<notin> X \<Longrightarrow> f x = g x"
   2.163 +  shows "integrable M f \<longleftrightarrow> integrable M g"
   2.164 +  unfolding integrable_iff_bounded
   2.165 +proof (rule conj_cong)
   2.166 +  { assume "f \<in> borel_measurable M" then have "g \<in> borel_measurable M"
   2.167 +      by (rule measurable_discrete_difference[where X=X]) (auto simp: assms) }
   2.168 +  moreover
   2.169 +  { assume "g \<in> borel_measurable M" then have "f \<in> borel_measurable M"
   2.170 +      by (rule measurable_discrete_difference[where X=X]) (auto simp: assms) }
   2.171 +  ultimately show "f \<in> borel_measurable M \<longleftrightarrow> g \<in> borel_measurable M" ..
   2.172 +next
   2.173 +  have "AE x in M. x \<notin> X"
   2.174 +    by (rule AE_discrete_difference) fact+
   2.175 +  then have "(\<integral>\<^sup>+ x. norm (f x) \<partial>M) = (\<integral>\<^sup>+ x. norm (g x) \<partial>M)"
   2.176 +    by (intro nn_integral_cong_AE) (auto simp: eq)
   2.177 +  then show "(\<integral>\<^sup>+ x. norm (f x) \<partial>M) < \<infinity> \<longleftrightarrow> (\<integral>\<^sup>+ x. norm (g x) \<partial>M) < \<infinity>"
   2.178 +    by simp
   2.179  qed
   2.180  
   2.181 +lemma integral_discrete_difference:
   2.182 +  fixes f :: "'a \<Rightarrow> 'b::{banach, second_countable_topology}"
   2.183 +  assumes X: "countable X"
   2.184 +  assumes null: "\<And>x. x \<in> X \<Longrightarrow> emeasure M {x} = 0" 
   2.185 +  assumes sets: "\<And>x. x \<in> X \<Longrightarrow> {x} \<in> sets M"
   2.186 +  assumes eq: "\<And>x. x \<in> space M \<Longrightarrow> x \<notin> X \<Longrightarrow> f x = g x"
   2.187 +  shows "integral\<^sup>L M f = integral\<^sup>L M g"
   2.188 +proof (rule integral_eq_cases)
   2.189 +  show eq: "integrable M f \<longleftrightarrow> integrable M g"
   2.190 +    by (rule integrable_discrete_difference[where X=X]) fact+
   2.191 +
   2.192 +  assume f: "integrable M f"
   2.193 +  show "integral\<^sup>L M f = integral\<^sup>L M g"
   2.194 +  proof (rule integral_cong_AE)
   2.195 +    show "f \<in> borel_measurable M" "g \<in> borel_measurable M"
   2.196 +      using f eq by (auto intro: borel_measurable_integrable)
   2.197 +
   2.198 +    have "AE x in M. x \<notin> X"
   2.199 +      by (rule AE_discrete_difference) fact+
   2.200 +    with AE_space show "AE x in M. f x = g x"
   2.201 +      by eventually_elim fact
   2.202 +  qed
   2.203 +qed
   2.204 +
   2.205 +lemma has_bochner_integral_discrete_difference:
   2.206 +  fixes f :: "'a \<Rightarrow> 'b::{banach, second_countable_topology}"
   2.207 +  assumes X: "countable X"
   2.208 +  assumes null: "\<And>x. x \<in> X \<Longrightarrow> emeasure M {x} = 0" 
   2.209 +  assumes sets: "\<And>x. x \<in> X \<Longrightarrow> {x} \<in> sets M"
   2.210 +  assumes eq: "\<And>x. x \<in> space M \<Longrightarrow> x \<notin> X \<Longrightarrow> f x = g x"
   2.211 +  shows "has_bochner_integral M f x \<longleftrightarrow> has_bochner_integral M g x"
   2.212 +  using integrable_discrete_difference[of X M f g, OF assms]
   2.213 +  using integral_discrete_difference[of X M f g, OF assms]
   2.214 +  by (metis has_bochner_integral_iff)
   2.215 +
   2.216 +lemma has_bochner_integral_even_function:
   2.217 +  fixes f :: "real \<Rightarrow> 'a :: {banach, second_countable_topology}"
   2.218 +  assumes f: "has_bochner_integral lborel (\<lambda>x. indicator {0..} x *\<^sub>R f x) x"
   2.219 +  assumes even: "\<And>x. f (- x) = f x"
   2.220 +  shows "has_bochner_integral lborel f (2 *\<^sub>R x)"
   2.221 +proof -
   2.222 +  have indicator: "\<And>x::real. indicator {..0} (- x) = indicator {0..} x"
   2.223 +    by (auto split: split_indicator)
   2.224 +  have "has_bochner_integral lborel (\<lambda>x. indicator {.. 0} x *\<^sub>R f x) x"
   2.225 +    by (subst lborel_has_bochner_integral_real_affine_iff[where c="-1" and t=0])
   2.226 +       (auto simp: indicator even f)
   2.227 +  with f have "has_bochner_integral lborel (\<lambda>x. indicator {0..} x *\<^sub>R f x + indicator {.. 0} x *\<^sub>R f x) (x + x)"
   2.228 +    by (rule has_bochner_integral_add)
   2.229 +  then have "has_bochner_integral lborel f (x + x)"
   2.230 +    by (rule has_bochner_integral_discrete_difference[where X="{0}", THEN iffD1, rotated 4])
   2.231 +       (auto split: split_indicator)
   2.232 +  then show ?thesis
   2.233 +    by (simp add: scaleR_2)
   2.234 +qed
   2.235 +
   2.236 +lemma has_bochner_integral_odd_function:
   2.237 +  fixes f :: "real \<Rightarrow> 'a :: {banach, second_countable_topology}"
   2.238 +  assumes f: "has_bochner_integral lborel (\<lambda>x. indicator {0..} x *\<^sub>R f x) x"
   2.239 +  assumes odd: "\<And>x. f (- x) = - f x"
   2.240 +  shows "has_bochner_integral lborel f 0"
   2.241 +proof -
   2.242 +  have indicator: "\<And>x::real. indicator {..0} (- x) = indicator {0..} x"
   2.243 +    by (auto split: split_indicator)
   2.244 +  have "has_bochner_integral lborel (\<lambda>x. - indicator {.. 0} x *\<^sub>R f x) x"
   2.245 +    by (subst lborel_has_bochner_integral_real_affine_iff[where c="-1" and t=0])
   2.246 +       (auto simp: indicator odd f)
   2.247 +  from has_bochner_integral_minus[OF this]
   2.248 +  have "has_bochner_integral lborel (\<lambda>x. indicator {.. 0} x *\<^sub>R f x) (- x)"
   2.249 +    by simp 
   2.250 +  with f have "has_bochner_integral lborel (\<lambda>x. indicator {0..} x *\<^sub>R f x + indicator {.. 0} x *\<^sub>R f x) (x + - x)"
   2.251 +    by (rule has_bochner_integral_add)
   2.252 +  then have "has_bochner_integral lborel f (x + - x)"
   2.253 +    by (rule has_bochner_integral_discrete_difference[where X="{0}", THEN iffD1, rotated 4])
   2.254 +       (auto split: split_indicator)
   2.255 +  then show ?thesis
   2.256 +    by simp
   2.257 +qed
   2.258 +
   2.259 +
   2.260  lemma filterlim_power2_at_top[intro]: "filterlim (\<lambda>(x::real). x\<^sup>2) at_top at_top"
   2.261    by (auto intro!: filterlim_at_top_mult_at_top filterlim_ident simp: power2_eq_square)
   2.262  
   2.263 @@ -897,6 +1086,7 @@
   2.264  
   2.265  subsection {* Normal distribution *}
   2.266  
   2.267 +
   2.268  definition normal_density :: "real \<Rightarrow> real \<Rightarrow> real \<Rightarrow> real" where
   2.269    "normal_density \<mu> \<sigma> x = 1 / sqrt (2 * pi * \<sigma>\<^sup>2) * exp (-(x - \<mu>)\<^sup>2/ (2 * \<sigma>\<^sup>2))"
   2.270  
   2.271 @@ -906,45 +1096,54 @@
   2.272  lemma std_normal_density_def: "std_normal_density x = (1 / sqrt (2 * pi)) * exp (- x\<^sup>2 / 2)"
   2.273    unfolding normal_density_def by simp
   2.274  
   2.275 +lemma normal_density_nonneg: "0 \<le> normal_density \<mu> \<sigma> x"
   2.276 +  by (auto simp: normal_density_def)
   2.277 +
   2.278 +lemma normal_density_pos: "0 < \<sigma> \<Longrightarrow> 0 < normal_density \<mu> \<sigma> x"
   2.279 +  by (auto simp: normal_density_def)
   2.280 +
   2.281  lemma borel_measurable_normal_density[measurable]: "normal_density \<mu> \<sigma> \<in> borel_measurable borel"
   2.282    by (auto simp: normal_density_def[abs_def])
   2.283  
   2.284 -lemma nn_integral_gaussian: "(\<integral>\<^sup>+x. (exp (- x\<^sup>2)) \<partial>lborel) = sqrt pi"
   2.285 +lemma gaussian_moment_0:
   2.286 +  "has_bochner_integral lborel (\<lambda>x. indicator {0..} x *\<^sub>R exp (- x\<^sup>2)) (sqrt pi / 2)"
   2.287  proof -
   2.288    let ?pI = "\<lambda>f. (\<integral>\<^sup>+s. f (s::real) * indicator {0..} s \<partial>lborel)"
   2.289    let ?gauss = "\<lambda>x. exp (- x\<^sup>2)"
   2.290  
   2.291 -  have "?pI ?gauss * ?pI ?gauss = ?pI (\<lambda>s. ?pI (\<lambda>x. ereal (x * exp (-x\<^sup>2 * (1 + s\<^sup>2)))))"
   2.292 -  proof-
   2.293 -    let ?ff= "\<lambda>(x, s). ((x * exp (- x\<^sup>2 * (1 + s\<^sup>2)) * indicator {0<..} s * indicator {0<..} x)) :: real"
   2.294 -  
   2.295 -    have *: "?pI ?gauss = (\<integral>\<^sup>+x. ?gauss x * indicator {0<..} x \<partial>lborel)"
   2.296 -      by (intro nn_integral_cong_AE AE_I[where N="{0}"]) (auto split: split_indicator)
   2.297 -  
   2.298 -    have "?pI ?gauss * ?pI ?gauss = (\<integral>\<^sup>+x. \<integral>\<^sup>+s. ?ff (x, s) \<partial>lborel \<partial>lborel)"
   2.299 -      unfolding *
   2.300 -      apply (auto simp: nn_integral_nonneg nn_integral_cmult[symmetric])
   2.301 -      apply (auto intro!: nn_integral_cong split:split_indicator)
   2.302 -      apply (auto simp: nn_integral_multc[symmetric])
   2.303 -      apply (subst nn_integral_real_affine[where t="0" and c="x"])
   2.304 -      by (auto simp: mult_exp_exp nn_integral_cmult[symmetric] field_simps zero_less_mult_iff
   2.305 -          intro!: nn_integral_cong split:split_indicator)
   2.306 -    also have "... = \<integral>\<^sup>+ (s::real). \<integral>\<^sup>+ (x::real). ?ff (x, s) \<partial>lborel \<partial>lborel"
   2.307 -      by (rule lborel_pair.Fubini[symmetric]) auto
   2.308 -    also have "... = ?pI (\<lambda>s. ?pI (\<lambda>x. ereal (x * exp (-x\<^sup>2 * (1 + s\<^sup>2)))))"
   2.309 -      by (rule nn_integral_cong_AE) (auto intro!: nn_integral_cong_AE AE_I[where N="{0}"] split: split_indicator_asm)
   2.310 -    finally show ?thesis
   2.311 -      by simp
   2.312 -  qed
   2.313 +  let ?I = "indicator {0<..} :: real \<Rightarrow> real"
   2.314 +  let ?ff= "\<lambda>x s. x * exp (- x\<^sup>2 * (1 + s\<^sup>2)) :: real"
   2.315 +
   2.316 +  have *: "?pI ?gauss = (\<integral>\<^sup>+x. ?gauss x * ?I x \<partial>lborel)"
   2.317 +    by (intro nn_integral_cong_AE AE_I[where N="{0}"]) (auto split: split_indicator)
   2.318 +
   2.319 +  have "?pI ?gauss * ?pI ?gauss = (\<integral>\<^sup>+x. \<integral>\<^sup>+s. ?gauss x * ?gauss s * ?I s * ?I x \<partial>lborel \<partial>lborel)"
   2.320 +    by (auto simp: nn_integral_nonneg nn_integral_cmult[symmetric] nn_integral_multc[symmetric] *
   2.321 +             intro!: nn_integral_cong split: split_indicator)
   2.322 +  also have "\<dots> = (\<integral>\<^sup>+x. \<integral>\<^sup>+s. ?ff x s * ?I s * ?I x \<partial>lborel \<partial>lborel)"
   2.323 +  proof (rule nn_integral_cong, cases)
   2.324 +    fix x :: real assume "x \<noteq> 0"
   2.325 +    then show "(\<integral>\<^sup>+s. ?gauss x * ?gauss s * ?I s * ?I x \<partial>lborel) = (\<integral>\<^sup>+s. ?ff x s * ?I s * ?I x \<partial>lborel)"
   2.326 +      by (subst nn_integral_real_affine[where t="0" and c="x"])
   2.327 +         (auto simp: mult_exp_exp nn_integral_cmult[symmetric] field_simps zero_less_mult_iff
   2.328 +               intro!: nn_integral_cong split: split_indicator)
   2.329 +  qed simp
   2.330 +  also have "... = \<integral>\<^sup>+s. \<integral>\<^sup>+x. ?ff x s * ?I s * ?I x \<partial>lborel \<partial>lborel"
   2.331 +    by (rule lborel_pair.Fubini'[symmetric]) auto
   2.332 +  also have "... = ?pI (\<lambda>s. ?pI (\<lambda>x. ?ff x s))"
   2.333 +    by (rule nn_integral_cong_AE)
   2.334 +       (auto intro!: nn_integral_cong_AE AE_I[where N="{0}"] split: split_indicator_asm)
   2.335    also have "\<dots> = ?pI (\<lambda>s. ereal (1 / (2 * (1 + s\<^sup>2))))"
   2.336    proof (intro nn_integral_cong ereal_right_mult_cong)
   2.337 -    fix s :: real show "?pI (\<lambda>x. ereal (x * exp (-x\<^sup>2 * (1 + s\<^sup>2)))) = ereal (1 / (2 * (1 + s\<^sup>2)))"
   2.338 +    fix s :: real show "?pI (\<lambda>x. ?ff x s) = ereal (1 / (2 * (1 + s\<^sup>2)))"
   2.339      proof (subst nn_integral_FTC_atLeast)
   2.340        have "((\<lambda>a. - (exp (- (a\<^sup>2 * (1 + s\<^sup>2))) / (2 + 2 * s\<^sup>2))) ---> (- (0 / (2 + 2 * s\<^sup>2)))) at_top"
   2.341          apply (intro tendsto_intros filterlim_compose[OF exp_at_bot] filterlim_compose[OF filterlim_uminus_at_bot_at_top])
   2.342          apply (subst mult_commute)         
   2.343 -        by (auto intro!: filterlim_tendsto_pos_mult_at_top filterlim_at_top_mult_at_top[OF filterlim_ident filterlim_ident] 
   2.344 -          simp: add_pos_nonneg  power2_eq_square add_nonneg_eq_0_iff)
   2.345 +        apply (auto intro!: filterlim_tendsto_pos_mult_at_top
   2.346 +                            filterlim_at_top_mult_at_top[OF filterlim_ident filterlim_ident] 
   2.347 +                    simp: add_pos_nonneg  power2_eq_square add_nonneg_eq_0_iff)
   2.348 +        done
   2.349        then show "((\<lambda>a. - (exp (- a\<^sup>2 - s\<^sup>2 * a\<^sup>2) / (2 + 2 * s\<^sup>2))) ---> 0) at_top"
   2.350          by (simp add: field_simps)
   2.351      qed (auto intro!: derivative_eq_intros simp: field_simps add_nonneg_eq_0_iff)
   2.352 @@ -958,47 +1157,307 @@
   2.353      by (simp add: power2_eq_square)
   2.354    then have "?pI ?gauss = sqrt (pi / 4)"
   2.355      using power_eq_iff_eq_base[of 2 "real (?pI ?gauss)" "sqrt (pi / 4)"]
   2.356 -      nn_integral_nonneg[of lborel "\<lambda>x. ?gauss x * indicator {0..} x"]
   2.357 +          nn_integral_nonneg[of lborel "\<lambda>x. ?gauss x * indicator {0..} x"]
   2.358      by (cases "?pI ?gauss") auto
   2.359 -  then show ?thesis
   2.360 -    by (subst nn_integral_even_function) (auto simp add: real_sqrt_divide)
   2.361 +  also have "?pI ?gauss = (\<integral>\<^sup>+x. indicator {0..} x *\<^sub>R exp (- x\<^sup>2) \<partial>lborel)"
   2.362 +    by (intro nn_integral_cong) (simp split: split_indicator)
   2.363 +  also have "sqrt (pi / 4) = sqrt pi / 2"
   2.364 +    by (simp add: real_sqrt_divide)
   2.365 +  finally show ?thesis
   2.366 +    by (rule has_bochner_integral_nn_integral[rotated 2]) auto
   2.367 +qed
   2.368 +
   2.369 +lemma gaussian_moment_1:
   2.370 +  "has_bochner_integral lborel (\<lambda>x::real. indicator {0..} x *\<^sub>R (exp (- x\<^sup>2) * x)) (1 / 2)" 
   2.371 +proof - 
   2.372 +  have "(\<integral>\<^sup>+x. indicator {0..} x *\<^sub>R (exp (- x\<^sup>2) * x) \<partial>lborel) =
   2.373 +    (\<integral>\<^sup>+x. ereal (x * exp (- x\<^sup>2)) * indicator {0..} x \<partial>lborel)"
   2.374 +    by (intro nn_integral_cong)
   2.375 +       (auto simp: ac_simps times_ereal.simps(1)[symmetric] ereal_indicator simp del: times_ereal.simps)
   2.376 +  also have "\<dots> = ereal (0 - (- exp (- 0\<^sup>2) / 2))"
   2.377 +  proof (rule nn_integral_FTC_atLeast)
   2.378 +    have "((\<lambda>x::real. - exp (- x\<^sup>2) / 2) ---> - 0 / 2) at_top"
   2.379 +      by (intro tendsto_divide tendsto_minus filterlim_compose[OF exp_at_bot]
   2.380 +                   filterlim_compose[OF filterlim_uminus_at_bot_at_top])
   2.381 +         auto
   2.382 +    then show "((\<lambda>a::real. - exp (- a\<^sup>2) / 2) ---> 0) at_top"
   2.383 +      by simp
   2.384 +  qed (auto intro!: derivative_eq_intros)
   2.385 +  also have "\<dots> = ereal (1 / 2)"
   2.386 +    by simp
   2.387 +  finally show ?thesis
   2.388 +    by (rule has_bochner_integral_nn_integral[rotated 2]) (auto split: split_indicator)
   2.389  qed
   2.390  
   2.391 -lemma has_bochner_integral_gaussian: "has_bochner_integral lborel (\<lambda>x. exp (- x\<^sup>2)) (sqrt pi)"
   2.392 -  by (auto intro!: nn_integral_gaussian has_bochner_integral_nn_integral)
   2.393 +lemma
   2.394 +  fixes k :: nat
   2.395 +  shows gaussian_moment_even_pos:
   2.396 +    "has_bochner_integral lborel (\<lambda>x. indicator {0..} x *\<^sub>R (exp (-x\<^sup>2)*x^(2 * k)))
   2.397 +       ((sqrt pi / 2) * (fact (2 * k) / (2 ^ (2 * k) * fact k)))" (is "?even")
   2.398 +    and gaussian_moment_odd_pos:
   2.399 +      "has_bochner_integral lborel (\<lambda>x. indicator {0..} x *\<^sub>R (exp (-x\<^sup>2)*x^(2 * k + 1))) (fact k / 2)" (is "?odd")
   2.400 +proof -
   2.401 +  let ?M = "\<lambda>k x. exp (- x\<^sup>2) * x^k :: real"
   2.402 +
   2.403 +  { fix k I assume Mk: "has_bochner_integral lborel (\<lambda>x. indicator {0..} x *\<^sub>R ?M k x) I"
   2.404 +    have "2 \<noteq> (0::real)"
   2.405 +      by linarith
   2.406 +    let ?f = "\<lambda>b. \<integral>x. indicator {0..} x *\<^sub>R ?M (k + 2) x * indicator {..b} x \<partial>lborel"
   2.407 +    have "((\<lambda>b. (k + 1) / 2 * (\<integral>x. indicator {..b} x *\<^sub>R (indicator {0..} x *\<^sub>R ?M k x) \<partial>lborel) - ?M (k + 1) b / 2) --->
   2.408 +        (k + 1) / 2 * (\<integral>x. indicator {0..} x *\<^sub>R ?M k x \<partial>lborel) - 0 / 2) at_top" (is ?tendsto)
   2.409 +    proof (intro tendsto_intros `2 \<noteq> 0` tendsto_integral_at_top sets_lborel Mk[THEN integrable.intros])
   2.410 +      show "(?M (k + 1) ---> 0) at_top"
   2.411 +      proof cases
   2.412 +        assume "even k"
   2.413 +        have "((\<lambda>x. ((x\<^sup>2)^(k div 2 + 1) / exp (x\<^sup>2)) * (1 / x) :: real) ---> 0 * 0) at_top"
   2.414 +          by (intro tendsto_intros tendsto_divide_0[OF tendsto_const] filterlim_compose[OF tendsto_power_div_exp_0]
   2.415 +                   filterlim_at_top_imp_at_infinity filterlim_ident filterlim_power2_at_top)
   2.416 +        also have "(\<lambda>x. ((x\<^sup>2)^(k div 2 + 1) / exp (x\<^sup>2)) * (1 / x) :: real) = ?M (k + 1)"
   2.417 +          using `even k` by (auto simp: even_mult_two_ex fun_eq_iff exp_minus field_simps power2_eq_square power_mult)
   2.418 +        finally show ?thesis by simp
   2.419 +      next
   2.420 +        assume "odd k"
   2.421 +        have "((\<lambda>x. ((x\<^sup>2)^((k - 1) div 2 + 1) / exp (x\<^sup>2)) :: real) ---> 0) at_top"
   2.422 +          by (intro filterlim_compose[OF tendsto_power_div_exp_0] filterlim_at_top_imp_at_infinity filterlim_ident filterlim_power2_at_top)
   2.423 +        also have "(\<lambda>x. ((x\<^sup>2)^((k - 1) div 2 + 1) / exp (x\<^sup>2)) :: real) = ?M (k + 1)"
   2.424 +          using `odd k` by (auto simp: odd_Suc_mult_two_ex fun_eq_iff exp_minus field_simps power2_eq_square power_mult)
   2.425 +        finally show ?thesis by simp
   2.426 +      qed
   2.427 +    qed
   2.428 +    also have "?tendsto \<longleftrightarrow> ((?f ---> (k + 1) / 2 * (\<integral>x. indicator {0..} x *\<^sub>R ?M k x \<partial>lborel) - 0 / 2) at_top)"
   2.429 +    proof (intro filterlim_cong refl eventually_at_top_linorder[THEN iffD2] exI[of _ 0] allI impI)
   2.430 +      fix b :: real assume b: "0 \<le> b"
   2.431 +      have "Suc k * (\<integral>x. indicator {0..b} x *\<^sub>R ?M k x \<partial>lborel) = (\<integral>x. indicator {0..b} x *\<^sub>R (exp (- x\<^sup>2) * ((Suc k) * x ^ k)) \<partial>lborel)"
   2.432 +        unfolding integral_mult_right_zero[symmetric] by (intro integral_cong) auto
   2.433 +      also have "\<dots> = exp (- b\<^sup>2) * b ^ (Suc k) - exp (- 0\<^sup>2) * 0 ^ (Suc k) -
   2.434 +          (\<integral>x. indicator {0..b} x *\<^sub>R (- 2 * x * exp (- x\<^sup>2) * x ^ (Suc k)) \<partial>lborel)"
   2.435 +        by (rule integral_by_parts')
   2.436 +           (auto intro!: derivative_eq_intros b
   2.437 +                 simp: real_of_nat_def[symmetric] diff_Suc real_of_nat_Suc field_simps split: nat.split)
   2.438 +      also have "(\<integral>x. indicator {0..b} x *\<^sub>R (- 2 * x * exp (- x\<^sup>2) * x ^ (Suc k)) \<partial>lborel) =
   2.439 +        (\<integral>x. indicator {0..b} x *\<^sub>R (- 2 * (exp (- x\<^sup>2) * x ^ (k + 2))) \<partial>lborel)"
   2.440 +        by (intro integral_cong) auto
   2.441 +      finally have "Suc k * (\<integral>x. indicator {0..b} x *\<^sub>R ?M k x \<partial>lborel) =
   2.442 +        exp (- b\<^sup>2) * b ^ (Suc k) + 2 * (\<integral>x. indicator {0..b} x *\<^sub>R ?M (k + 2) x \<partial>lborel)"
   2.443 +        apply (simp del: real_scaleR_def integral_mult_right add: integral_mult_right[symmetric])
   2.444 +        apply (subst integral_mult_right_zero[symmetric])
   2.445 +        apply (intro integral_cong)
   2.446 +        apply simp_all
   2.447 +        done
   2.448 +      then show "(k + 1) / 2 * (\<integral>x. indicator {..b} x *\<^sub>R (indicator {0..} x *\<^sub>R ?M k x)\<partial>lborel) - exp (- b\<^sup>2) * b ^ (k + 1) / 2 = ?f b"
   2.449 +        by (simp add: field_simps atLeastAtMost_def indicator_inter_arith)
   2.450 +    qed
   2.451 +    finally have int_M_at_top: "((?f ---> (k + 1) / 2 * (\<integral>x. indicator {0..} x *\<^sub>R ?M k x \<partial>lborel)) at_top)"
   2.452 +      by simp
   2.453 +    
   2.454 +    have "has_bochner_integral lborel (\<lambda>x. indicator {0..} x *\<^sub>R ?M (k + 2) x) ((k + 1) / 2 * I)"
   2.455 +    proof (rule has_bochner_integral_monotone_convergence_at_top)
   2.456 +      fix y :: real
   2.457 +      have *: "(\<lambda>x. indicator {0..} x *\<^sub>R ?M (k + 2) x * indicator {..y} x::real) =
   2.458 +            (\<lambda>x. indicator {0..y} x *\<^sub>R ?M (k + 2) x)"
   2.459 +        by rule (simp split: split_indicator)
   2.460 +      show "integrable lborel (\<lambda>x. indicator {0..} x *\<^sub>R (?M (k + 2) x) * indicator {..y} x::real)"
   2.461 +        unfolding * by (rule borel_integrable_compact) (auto intro!: continuous_intros)
   2.462 +      show "((?f ---> (k + 1) / 2 * I) at_top)"
   2.463 +        using int_M_at_top has_bochner_integral_integral_eq[OF Mk] by simp
   2.464 +    qed (auto split: split_indicator) }
   2.465 +  note step = this
   2.466 +
   2.467 +  show ?even
   2.468 +  proof (induct k)
   2.469 +    case (Suc k)
   2.470 +    note step[OF this]
   2.471 +    also have "(real (2 * k + 1) / 2 * (sqrt pi / 2 * (real (fact (2 * k)) / real (2 ^ (2 * k) * fact k)))) =
   2.472 +      sqrt pi / 2 * (real (fact (2 * Suc k)) / real (2 ^ (2 * Suc k) * fact (Suc k)))"
   2.473 +      by (simp add: field_simps real_of_nat_Suc divide_simps del: fact_Suc) (simp add: field_simps)
   2.474 +    finally show ?case
   2.475 +      by simp
   2.476 +  qed (insert gaussian_moment_0, simp)
   2.477  
   2.478 -lemma integral_gaussian: "(\<integral>x. (exp (- x\<^sup>2)) \<partial>lborel) = sqrt pi"
   2.479 -  using has_bochner_integral_gaussian by (rule has_bochner_integral_integral_eq)
   2.480 +  show ?odd
   2.481 +  proof (induct k)
   2.482 +    case (Suc k)
   2.483 +    note step[OF this]
   2.484 +    also have "(real (2 * k + 1 + 1) / 2 * (real (fact k) / 2)) = real (fact (Suc k)) / 2"
   2.485 +      by (simp add: field_simps real_of_nat_Suc divide_simps del: fact_Suc) (simp add: field_simps)
   2.486 +    finally show ?case
   2.487 +      by simp
   2.488 +  qed (insert gaussian_moment_1, simp)
   2.489 +qed
   2.490 +
   2.491 +context
   2.492 +  fixes k :: nat and \<mu> \<sigma> :: real assumes [arith]: "0 < \<sigma>"
   2.493 +begin
   2.494 +
   2.495 +lemma normal_moment_even:
   2.496 +  "has_bochner_integral lborel (\<lambda>x. normal_density \<mu> \<sigma> x * (x - \<mu>) ^ (2 * k)) (fact (2 * k) / ((2 / \<sigma>\<^sup>2)^k * fact k))"
   2.497 +proof -
   2.498 +  have eq: "\<And>x::real. x\<^sup>2^k = (x^k)\<^sup>2"
   2.499 +    by (simp add: power_mult[symmetric] ac_simps)
   2.500 +
   2.501 +  have "has_bochner_integral lborel (\<lambda>x. exp (-x\<^sup>2)*x^(2 * k))
   2.502 +      (sqrt pi * (fact (2 * k) / (2 ^ (2 * k) * fact k)))"
   2.503 +    using has_bochner_integral_even_function[OF gaussian_moment_even_pos[where k=k]] by simp
   2.504 +  then have "has_bochner_integral lborel (\<lambda>x. (exp (-x\<^sup>2)*x^(2 * k)) * ((2*\<sigma>\<^sup>2)^k / sqrt (2 * pi * \<sigma>\<^sup>2)))
   2.505 +      ((sqrt pi * (fact (2 * k) / (2 ^ (2 * k) * fact k))) * ((2*\<sigma>\<^sup>2)^k / sqrt (2 * pi * \<sigma>\<^sup>2)))"
   2.506 +    by (rule has_bochner_integral_mult_left)
   2.507 +  also have "(\<lambda>x. (exp (-x\<^sup>2)*x^(2 * k)) * ((2*\<sigma>\<^sup>2)^k / sqrt (2 * pi * \<sigma>\<^sup>2))) =
   2.508 +    (\<lambda>x. exp (- ((sqrt 2 * \<sigma>) * x)\<^sup>2 / (2*\<sigma>\<^sup>2)) * ((sqrt 2 * \<sigma>) * x) ^ (2 * k) / sqrt (2 * pi * \<sigma>\<^sup>2))"
   2.509 +    by (auto simp: fun_eq_iff field_simps real_sqrt_power[symmetric] real_sqrt_mult
   2.510 +                   real_sqrt_divide power_mult eq)
   2.511 +  also have "((sqrt pi * (fact (2 * k) / (2 ^ (2 * k) * fact k))) * ((2*\<sigma>\<^sup>2)^k / sqrt (2 * pi * \<sigma>\<^sup>2))) = 
   2.512 +    (inverse (sqrt 2 * \<sigma>) * (real (fact (2 * k))) / ((2/\<sigma>\<^sup>2) ^ k * real (fact k)))"
   2.513 +    by (auto simp: fun_eq_iff power_mult field_simps real_sqrt_power[symmetric] real_sqrt_mult
   2.514 +                   power2_eq_square)
   2.515 +  finally show ?thesis
   2.516 +    unfolding normal_density_def
   2.517 +    by (subst lborel_has_bochner_integral_real_affine_iff[where c="sqrt 2 * \<sigma>" and t=\<mu>]) simp_all
   2.518 +qed
   2.519  
   2.520 -lemma integrable_gaussian[intro]: "integrable lborel (\<lambda>x. exp (- x\<^sup>2)::real)"
   2.521 -  using has_bochner_integral_gaussian by rule
   2.522 +lemma normal_moment_abs_odd:
   2.523 +  "has_bochner_integral lborel (\<lambda>x. normal_density \<mu> \<sigma> x * \<bar>x - \<mu>\<bar>^(2 * k + 1)) (2^k * \<sigma>^(2 * k + 1) * fact k * sqrt (2 / pi))"
   2.524 +proof -
   2.525 +  have "has_bochner_integral lborel (\<lambda>x. indicator {0..} x *\<^sub>R (exp (-x\<^sup>2)*\<bar>x\<bar>^(2 * k + 1))) (fact k / 2)"
   2.526 +    by (rule has_bochner_integral_cong[THEN iffD1, OF _ _ _ gaussian_moment_odd_pos[of k]]) auto
   2.527 +  from has_bochner_integral_even_function[OF this]
   2.528 +  have "has_bochner_integral lborel (\<lambda>x. exp (-x\<^sup>2)*\<bar>x\<bar>^(2 * k + 1)) (fact k)"
   2.529 +    by simp
   2.530 +  then have "has_bochner_integral lborel (\<lambda>x. (exp (-x\<^sup>2)*\<bar>x\<bar>^(2 * k + 1)) * (2^k * \<sigma>^(2 * k + 1) / sqrt (pi * \<sigma>\<^sup>2)))
   2.531 +      (fact k * (2^k * \<sigma>^(2 * k + 1) / sqrt (pi * \<sigma>\<^sup>2)))"
   2.532 +    by (rule has_bochner_integral_mult_left)
   2.533 +  also have "(\<lambda>x. (exp (-x\<^sup>2)*\<bar>x\<bar>^(2 * k + 1)) * (2^k * \<sigma>^(2 * k + 1) / sqrt (pi * \<sigma>\<^sup>2))) =
   2.534 +    (\<lambda>x. exp (- (((sqrt 2 * \<sigma>) * x)\<^sup>2 / (2 * \<sigma>\<^sup>2))) * \<bar>sqrt 2 * \<sigma> * x\<bar> ^ (2 * k + 1) / sqrt (2 * pi * \<sigma>\<^sup>2))"
   2.535 +    by (simp add: field_simps abs_mult real_sqrt_power[symmetric] power_mult real_sqrt_mult)
   2.536 +  also have "(fact k * (2^k * \<sigma>^(2 * k + 1) / sqrt (pi * \<sigma>\<^sup>2))) = 
   2.537 +    (inverse (sqrt 2) * inverse \<sigma> * (2 ^ k * (\<sigma> * \<sigma> ^ (2 * k)) * real (fact k) * sqrt (2 / pi)))"
   2.538 +    by (auto simp: fun_eq_iff power_mult field_simps real_sqrt_power[symmetric] real_sqrt_divide
   2.539 +                   real_sqrt_mult)
   2.540 +  finally show ?thesis
   2.541 +    unfolding normal_density_def
   2.542 +    by (subst lborel_has_bochner_integral_real_affine_iff[where c="sqrt 2 * \<sigma>" and t=\<mu>])
   2.543 +       simp_all
   2.544 +qed
   2.545 +
   2.546 +lemma normal_moment_odd:
   2.547 +  "has_bochner_integral lborel (\<lambda>x. normal_density \<mu> \<sigma> x * (x - \<mu>)^(2 * k + 1)) 0"
   2.548 +proof -
   2.549 +  have "has_bochner_integral lborel (\<lambda>x. exp (- x\<^sup>2) * x^(2 * k + 1)::real) 0"
   2.550 +    using gaussian_moment_odd_pos by (rule has_bochner_integral_odd_function) simp
   2.551 +  then have "has_bochner_integral lborel (\<lambda>x. (exp (-x\<^sup>2)*x^(2 * k + 1)) * (2^k*\<sigma>^(2*k)/sqrt pi))
   2.552 +      (0 * (2^k*\<sigma>^(2*k)/sqrt pi))"
   2.553 +    by (rule has_bochner_integral_mult_left)
   2.554 +  also have "(\<lambda>x. (exp (-x\<^sup>2)*x^(2 * k + 1)) * (2^k*\<sigma>^(2*k)/sqrt pi)) =
   2.555 +    (\<lambda>x. exp (- ((sqrt 2 * \<sigma> * x)\<^sup>2 / (2 * \<sigma>\<^sup>2))) *
   2.556 +          (sqrt 2 * \<sigma> * x * (sqrt 2 * \<sigma> * x) ^ (2 * k)) /
   2.557 +          sqrt (2 * pi * \<sigma>\<^sup>2))"
   2.558 +    unfolding real_sqrt_mult
   2.559 +    by (simp add: field_simps abs_mult real_sqrt_power[symmetric] power_mult fun_eq_iff)
   2.560 +  finally show ?thesis
   2.561 +    unfolding normal_density_def
   2.562 +    by (subst lborel_has_bochner_integral_real_affine_iff[where c="sqrt 2 * \<sigma>" and t=\<mu>]) simp_all
   2.563 +qed
   2.564 +
   2.565 +lemma integral_normal_moment_even:
   2.566 +  "integral\<^sup>L lborel (\<lambda>x. normal_density \<mu> \<sigma> x * (x - \<mu>)^(2 * k)) = fact (2 * k) / ((2 / \<sigma>\<^sup>2)^k * fact k)"
   2.567 +  using normal_moment_even by (rule has_bochner_integral_integral_eq)
   2.568 +
   2.569 +lemma integral_normal_moment_abs_odd:
   2.570 +  "integral\<^sup>L lborel (\<lambda>x. normal_density \<mu> \<sigma> x * \<bar>x - \<mu>\<bar>^(2 * k + 1)) = 2 ^ k * \<sigma> ^ (2 * k + 1) * fact k * sqrt (2 / pi)"
   2.571 +  using normal_moment_abs_odd by (rule has_bochner_integral_integral_eq)
   2.572 +
   2.573 +lemma integral_normal_moment_odd:
   2.574 +  "integral\<^sup>L lborel (\<lambda>x. normal_density \<mu> \<sigma> x * (x - \<mu>)^(2 * k + 1)) = 0"
   2.575 +  using normal_moment_odd by (rule has_bochner_integral_integral_eq)
   2.576 +
   2.577 +end
   2.578 +
   2.579  
   2.580  context
   2.581    fixes \<sigma> :: real
   2.582    assumes \<sigma>_pos[arith]: "0 < \<sigma>"
   2.583  begin
   2.584  
   2.585 -lemma nn_integral_normal_density: "(\<integral>\<^sup>+x. normal_density \<mu> \<sigma> x \<partial>lborel) = 1"
   2.586 -  unfolding normal_density_def
   2.587 -  apply (subst times_ereal.simps(1)[symmetric],subst nn_integral_cmult)
   2.588 -  apply auto
   2.589 -  apply (subst nn_integral_real_affine[where t=\<mu> and  c="(sqrt 2) * \<sigma>"])
   2.590 -  by (auto simp: power_mult_distrib nn_integral_gaussian real_sqrt_mult one_ereal_def)
   2.591 +lemma normal_moment_nz_1: "has_bochner_integral lborel (\<lambda>x. normal_density \<mu> \<sigma> x * x) \<mu>"
   2.592 +proof -
   2.593 +  note normal_moment_even[OF \<sigma>_pos, of \<mu> 0]
   2.594 +  note normal_moment_odd[OF \<sigma>_pos, of \<mu> 0] has_bochner_integral_mult_left[of \<mu>, OF this]
   2.595 +  note has_bochner_integral_add[OF this]
   2.596 +  then show ?thesis
   2.597 +    by (simp add: power2_eq_square field_simps)  
   2.598 +qed
   2.599 +
   2.600 +lemma integral_normal_moment_nz_1:
   2.601 +  "integral\<^sup>L lborel (\<lambda>x. normal_density \<mu> \<sigma> x * x) = \<mu>"
   2.602 +  using normal_moment_nz_1 by (rule has_bochner_integral_integral_eq)
   2.603 +
   2.604 +lemma integrable_normal_moment_nz_1: "integrable lborel (\<lambda>x. normal_density \<mu> \<sigma> x * x)"
   2.605 +  using normal_moment_nz_1 by rule
   2.606  
   2.607 -lemma 
   2.608 -  shows normal_density_pos: "\<And>x. 0 < normal_density \<mu> \<sigma> x"
   2.609 -  and normal_density_nonneg: "\<And>x. 0 \<le> normal_density \<mu> \<sigma> x" 
   2.610 -  by (auto simp: normal_density_def)
   2.611 +lemma integrable_normal_moment: "integrable lborel (\<lambda>x. normal_density \<mu> \<sigma> x * (x - \<mu>)^k)"
   2.612 +proof cases
   2.613 +  assume "even k" then show ?thesis
   2.614 +    using integrable.intros[OF normal_moment_even] by (auto simp add: even_mult_two_ex)
   2.615 +next
   2.616 +  assume "odd k" then show ?thesis
   2.617 +    using integrable.intros[OF normal_moment_odd] by (auto simp add: odd_Suc_mult_two_ex)
   2.618 +qed
   2.619  
   2.620 -lemma integrable_normal[intro]: "integrable lborel (normal_density \<mu> \<sigma>)"
   2.621 -  by (auto intro!: integrableI_nn_integral_finite simp: nn_integral_normal_density normal_density_nonneg)
   2.622 +lemma integrable_normal_moment_abs: "integrable lborel (\<lambda>x. normal_density \<mu> \<sigma> x * \<bar>x - \<mu>\<bar>^k)"
   2.623 +proof cases
   2.624 +  assume "even k" then show ?thesis
   2.625 +    using integrable.intros[OF normal_moment_even] by (auto simp add: even_mult_two_ex power_even_abs)
   2.626 +next
   2.627 +  assume "odd k" then show ?thesis
   2.628 +    using integrable.intros[OF normal_moment_abs_odd] by (auto simp add: odd_Suc_mult_two_ex)
   2.629 +qed
   2.630 +
   2.631 +lemma integrable_normal_density[simp, intro]: "integrable lborel (normal_density \<mu> \<sigma>)"
   2.632 +  using integrable_normal_moment[of \<mu> 0] by simp
   2.633  
   2.634  lemma integral_normal_density[simp]: "(\<integral>x. normal_density \<mu> \<sigma> x \<partial>lborel) = 1"
   2.635 -  by (simp add: integral_eq_nn_integral normal_density_nonneg nn_integral_normal_density)
   2.636 +  using integral_normal_moment_even[of \<sigma> \<mu> 0] by simp
   2.637  
   2.638  lemma prob_space_normal_density:
   2.639 -  "prob_space (density lborel (normal_density \<mu> \<sigma>))" (is "prob_space ?D")
   2.640 -  proof qed (simp add: emeasure_density nn_integral_normal_density)
   2.641 +  "prob_space (density lborel (normal_density \<mu> \<sigma>))"
   2.642 +  proof qed (simp add: emeasure_density nn_integral_eq_integral normal_density_nonneg)
   2.643 +
   2.644 +end
   2.645 +
   2.646 +
   2.647 +
   2.648 +context
   2.649 +  fixes k :: nat
   2.650 +begin
   2.651 +
   2.652 +lemma std_normal_moment_even:
   2.653 +  "has_bochner_integral lborel (\<lambda>x. std_normal_density x * x ^ (2 * k)) (fact (2 * k) / (2^k * fact k))"
   2.654 +  using normal_moment_even[of 1 0 k] by simp
   2.655 +
   2.656 +lemma std_normal_moment_abs_odd:
   2.657 +  "has_bochner_integral lborel (\<lambda>x. std_normal_density x * \<bar>x\<bar>^(2 * k + 1)) (sqrt (2/pi) * 2^k * fact k)"
   2.658 +  using normal_moment_abs_odd[of 1 0 k] by (simp add: ac_simps)
   2.659 +
   2.660 +lemma std_normal_moment_odd:
   2.661 +  "has_bochner_integral lborel (\<lambda>x. std_normal_density x * x^(2 * k + 1)) 0"
   2.662 +  using normal_moment_odd[of 1 0 k] by simp
   2.663 +
   2.664 +lemma integral_std_normal_moment_even:
   2.665 +  "integral\<^sup>L lborel (\<lambda>x. std_normal_density x * x^(2*k)) = fact (2 * k) / (2^k * fact k)"
   2.666 +  using std_normal_moment_even by (rule has_bochner_integral_integral_eq)
   2.667 +
   2.668 +lemma integral_std_normal_moment_abs_odd:
   2.669 +  "integral\<^sup>L lborel (\<lambda>x. std_normal_density x * \<bar>x\<bar>^(2 * k + 1)) = sqrt (2 / pi) * 2^k * fact k"
   2.670 +  using std_normal_moment_abs_odd by (rule has_bochner_integral_integral_eq)
   2.671 +
   2.672 +lemma integral_std_normal_moment_odd:
   2.673 +  "integral\<^sup>L lborel (\<lambda>x. std_normal_density x * x^(2 * k + 1)) = 0"
   2.674 +  using std_normal_moment_odd by (rule has_bochner_integral_integral_eq)
   2.675 +
   2.676 +lemma integrable_std_normal_moment_abs: "integrable lborel (\<lambda>x. std_normal_density x * \<bar>x\<bar>^k)"
   2.677 +  using integrable_normal_moment_abs[of 1 0 k] by simp
   2.678 +
   2.679 +lemma integrable_std_normal_moment: "integrable lborel (\<lambda>x. std_normal_density x * x^k)"
   2.680 +  using integrable_normal_moment[of 1 0 k] by simp
   2.681  
   2.682  end
   2.683  
   2.684 @@ -1058,7 +1517,7 @@
   2.685      by (subst nn_integral_cmult[symmetric]) (auto simp: \<sigma>'_def \<tau>'_def normal_density_def)
   2.686  
   2.687    also have "... = (\<lambda>x. (normal_density 0 (sqrt (\<sigma>\<^sup>2 + \<tau>\<^sup>2)) x))"
   2.688 -    by (subst nn_integral_normal_density) (auto simp: sum_power2_gt_zero_iff)
   2.689 +    by (subst nn_integral_eq_integral) (auto simp: normal_density_nonneg)
   2.690  
   2.691    finally show ?thesis by fast
   2.692  qed
   2.693 @@ -1145,291 +1604,31 @@
   2.694      show ?case using 1 2 3 by simp  
   2.695  qed
   2.696  
   2.697 -lemma nn_integral_x_exp_x_square: "(\<integral>\<^sup>+x. ereal (x * exp (- x\<^sup>2 )) \<partial>lborel) = ereal 1 / 2" 
   2.698 -  and nn_integral_x_exp_x_square_indicator: "(\<integral>\<^sup>+x. ereal( x * exp (-x\<^sup>2 )) * indicator {0..} x \<partial>lborel) = ereal 1 / 2"
   2.699 -proof - 
   2.700 -  let ?F = "\<lambda>x. - exp (-x\<^sup>2 ) / 2"
   2.701 -
   2.702 -  have 1: "(\<integral>\<^sup>+x. ereal (x * exp (- x\<^sup>2)) * indicator {0..} x \<partial>lborel) =ereal( 0 - ?F 0)"
   2.703 -  apply (rule nn_integral_FTC_atLeast)
   2.704 -  apply (auto intro!: derivative_eq_intros)
   2.705 -  apply (rule tendsto_minus_cancel)
   2.706 -  apply (simp add: field_simps)
   2.707 -  proof - 
   2.708 -    have "((\<lambda>(x::real). exp (- x\<^sup>2) / 2) ---> 0 / 2) at_top"
   2.709 -    apply (intro tendsto_divide filterlim_compose[OF exp_at_bot] filterlim_compose[OF filterlim_uminus_at_bot_at_top])
   2.710 -    apply auto
   2.711 -    done
   2.712 -    then show "((\<lambda>(x::real). exp (- x\<^sup>2) / 2) ---> 0) at_top" by simp
   2.713 -  qed
   2.714 -
   2.715 -  also have 2: "(\<integral>\<^sup>+x. ereal (x * exp (- x\<^sup>2)) * indicator {0..} x \<partial>lborel) = integral\<^sup>N lborel (\<lambda>x. ereal (x * exp (- x\<^sup>2)))"
   2.716 -    by (subst(2) nn_integral_max_0[symmetric])
   2.717 -       (auto intro!: nn_integral_cong split: split_indicator simp: max_def zero_le_mult_iff)
   2.718 -  finally show "(\<integral>\<^sup>+x. ereal (x * exp (- x\<^sup>2)) \<partial>lborel) = ereal 1 / 2" by auto
   2.719 -
   2.720 -  show "(\<integral>\<^sup>+x. ereal (x * exp (- x\<^sup>2)) * indicator {0..} x \<partial>lborel) = ereal 1 / 2" using 1 by auto
   2.721 -qed
   2.722 -
   2.723 -lemma borel_integral_x_times_standard_normal[intro]: "(\<integral>x. std_normal_density x * x \<partial>lborel) = 0" 
   2.724 -  and borel_integrable_x_times_standard_normal[intro]: "integrable lborel (\<lambda>x. std_normal_density x * x)"
   2.725 -  and borel_integral_x_times_standard_normal'[intro]: "(\<integral>x. x * std_normal_density x \<partial>lborel) = 0" 
   2.726 -  and borel_integrable_x_times_standard_normal'[intro]: "integrable lborel (\<lambda>x. x * std_normal_density x)"
   2.727 -proof -    
   2.728 -  have 0: "(\<integral>\<^sup>+x. ereal (x * std_normal_density x) \<partial>lborel) = \<integral>\<^sup>+x. ereal (x * std_normal_density x) * indicator {0..} x \<partial>lborel"
   2.729 -    apply (subst nn_integral_max_0[symmetric]) 
   2.730 -    unfolding max_def std_normal_density_def
   2.731 -    apply (auto intro!: nn_integral_cong split:split_indicator simp: zero_le_divide_iff zero_le_mult_iff)
   2.732 -    apply (metis not_le pi_gt_zero)
   2.733 -    done
   2.734 -
   2.735 -  have 1: "(\<integral>\<^sup>+x. ereal (- (x * std_normal_density x)) \<partial>lborel) = \<integral>\<^sup>+x. ereal (x * std_normal_density x) * indicator {0..} x \<partial>lborel"
   2.736 -    apply (subst(2) nn_integral_real_affine[where c = "-1" and t = 0])
   2.737 -    apply(auto simp: std_normal_density_def split: split_indicator)
   2.738 -    apply(subst nn_integral_max_0[symmetric]) 
   2.739 -    unfolding max_def std_normal_density_def
   2.740 -    apply (auto intro!: nn_integral_cong split: split_indicator
   2.741 -                simp: divide_le_0_iff mult_le_0_iff one_ereal_def[symmetric])
   2.742 -    apply (metis not_le pi_gt_zero)
   2.743 -    done
   2.744 -
   2.745 -  have 2: "sqrt pi / sqrt 2 * (\<integral>\<^sup>+x. ereal (x * std_normal_density x) * indicator {0..} x \<partial>lborel) = integral\<^sup>N lborel (\<lambda>x. ereal (x * exp (- x\<^sup>2)))"
   2.746 -    unfolding std_normal_density_def
   2.747 -    apply (subst nn_integral_real_affine[where c = "sqrt 2" and t = 0])
   2.748 -    apply (auto simp: power_mult_distrib split: split_indicator)
   2.749 -    apply (subst mult_assoc[symmetric])
   2.750 -    apply (subst nn_integral_cmult[symmetric])
   2.751 -    apply auto
   2.752 -    apply (subst(2) nn_integral_max_0[symmetric])
   2.753 -    apply (auto intro!: nn_integral_cong split: split_indicator simp: max_def zero_le_mult_iff real_sqrt_mult)
   2.754 -    done
   2.755 -
   2.756 -  have *: "(\<integral>\<^sup>+x. ereal (x * std_normal_density x) * indicator {0..} x \<partial>lborel) = sqrt 2 / sqrt pi *(integral\<^sup>N lborel (\<lambda>x. ereal (x * exp (- x\<^sup>2))))"
   2.757 -    apply (subst 2[symmetric])
   2.758 -    apply (subst mult_assoc[symmetric])
   2.759 -    apply (auto simp: field_simps one_ereal_def[symmetric])
   2.760 -    done
   2.761 -    
   2.762 -  show "(\<integral> x. x * std_normal_density x \<partial>lborel) = 0" "integrable lborel (\<lambda>x. x * std_normal_density x)"
   2.763 -    by (subst real_lebesgue_integral_def)
   2.764 -       (auto simp: 0 1 * nn_integral_x_exp_x_square real_integrable_def)
   2.765 -
   2.766 -  then show "(\<integral> x. std_normal_density x * x \<partial>lborel) = 0" "integrable lborel (\<lambda>x. std_normal_density x * x)"
   2.767 -    by (simp_all add:mult_commute)
   2.768 -qed
   2.769 -
   2.770  lemma (in prob_space) standard_normal_distributed_expectation:
   2.771 -  assumes D: "distributed M lborel X std_normal_density "
   2.772 +  assumes D: "distributed M lborel X std_normal_density"
   2.773    shows "expectation X = 0"
   2.774 +  using integral_std_normal_moment_odd[of 0]
   2.775    by (auto simp: distributed_integral[OF D, of "\<lambda>x. x", symmetric])
   2.776  
   2.777  lemma (in prob_space) normal_distributed_expectation:
   2.778 -  assumes pos_var[arith]: "0 < \<sigma>"
   2.779 +  assumes \<sigma>[arith]: "0 < \<sigma>"
   2.780    assumes D: "distributed M lborel X (normal_density \<mu> \<sigma>)"
   2.781    shows "expectation X = \<mu>"
   2.782 -proof -
   2.783 -  def X' \<equiv> "\<lambda>x. (X x - \<mu>) / \<sigma>"
   2.784 -  then have D1: "distributed M lborel X' std_normal_density"
   2.785 -    by (simp add: normal_standard_normal_convert[OF pos_var, of X \<mu>, symmetric] D)
   2.786 -  then have [simp]: "integrable M X'"
   2.787 -    by (rule distributed_integrable_var) auto
   2.788 -
   2.789 -  have "expectation X = expectation (\<lambda>x. \<mu> + \<sigma> * X' x)"
   2.790 -    by (simp add: X'_def)
   2.791 -  then show ?thesis
   2.792 -    by (simp add: prob_space standard_normal_distributed_expectation[OF D1])
   2.793 -qed
   2.794 -
   2.795 -lemma integral_xsquare_exp_xsquare: "(\<integral> x. (x\<^sup>2 * exp (-x\<^sup>2 )) \<partial>lborel) =  sqrt pi / 2"
   2.796 -  and integrable_xsquare_exp_xsquare: "integrable lborel (\<lambda>x. x\<^sup>2 * exp (- x\<^sup>2)::real)"
   2.797 -proof- 
   2.798 -  note filterlim_compose[OF exp_at_top, intro] filterlim_ident[intro]
   2.799 -
   2.800 -  let ?f = "(\<lambda>x. x * - exp (- x\<^sup>2) / 2 - 0 * - exp (- 0\<^sup>2) / 2 -
   2.801 -                 \<integral> xa. 1 * (- exp (- xa\<^sup>2) / 2) * indicator {0..x} xa \<partial>lborel)::real\<Rightarrow>real"
   2.802 -  let ?IFunc = "(\<lambda>z. \<integral>x. (x\<^sup>2 * exp (- x\<^sup>2)) * indicator {0 .. z} x \<partial>lborel)::real\<Rightarrow>real"
   2.803 -
   2.804 -
   2.805 -  from nn_integral_gaussian  
   2.806 -  have 1: "(\<integral>\<^sup>+xa. ereal (exp (- xa\<^sup>2)) * indicator {0..} xa \<partial>lborel) = ereal (sqrt pi) / ereal 2"
   2.807 -    apply (subst (asm) nn_integral_even_function)
   2.808 -    apply simp
   2.809 -    apply simp
   2.810 -    apply (cases "\<integral>\<^sup>+ xa. ereal (exp (- xa\<^sup>2)) * indicator {0..} xa \<partial>lborel")
   2.811 -    apply auto
   2.812 -    done
   2.813 -
   2.814 -  then have I: "(\<integral>xa. exp (- xa\<^sup>2) * indicator {0..} xa \<partial>lborel) = sqrt pi / 2"
   2.815 -    by (subst integral_eq_nn_integral) (auto simp: ereal_mult_indicator)
   2.816 -  
   2.817 -  have byparts: "?IFunc = (\<lambda>z. (if z < 0 then 0 else ?f z))"
   2.818 -    proof (intro HOL.ext, subst split_ifs, safe)
   2.819 -      fix z::real assume [arith]:" \<not> z < 0 "
   2.820 -  
   2.821 -      have "?IFunc z =  \<integral>x. (x * (x * exp (- x\<^sup>2))) * indicator {0 .. z} x \<partial>lborel"
   2.822 -        by(auto intro!: integral_cong split: split_indicator simp: power2_eq_square)
   2.823 -  
   2.824 -      also have "... = (\<lambda>x. x) z * (\<lambda>x. - exp (- x\<^sup>2 ) / 2) z - (\<lambda>x. x) 0 * (\<lambda>x. - exp (- x\<^sup>2) / 2) 0 
   2.825 -          -  \<integral>x. 1 * ( - exp (- x\<^sup>2) / 2) * indicator {0 .. z} x \<partial>lborel"
   2.826 -        by(rule integral_by_parts, auto intro!: derivative_eq_intros)  
   2.827 -      finally have "?IFunc z = ..." .
   2.828 -
   2.829 -      then show "?IFunc z = ?f z" by simp
   2.830 -    qed simp    
   2.831 -
   2.832 -  have [simp]: "(\<lambda>y. \<integral> x. x\<^sup>2 * exp (- x\<^sup>2) * indicator {0..} x * indicator {..y} x \<partial>lborel) = ?IFunc"
   2.833 -    by(auto intro!: integral_cong split:split_indicator)
   2.834 -
   2.835 -  have [intro]: "((\<lambda>(x::real). x * exp (- x\<^sup>2) / 2) ---> 0) at_top"
   2.836 -  proof -
   2.837 -    have "((\<lambda>(x::real). x * exp (- x\<^sup>2) / 2) ---> 0 / 2) at_top"
   2.838 -      apply (intro tendsto_divide filterlim_compose[OF exp_at_bot] filterlim_compose[OF filterlim_uminus_at_bot_at_top])
   2.839 -      apply (auto simp: exp_minus inverse_eq_divide)
   2.840 -      apply (rule lhospital_at_top_at_top[where f' = "\<lambda>x. 1" and g' = "\<lambda>x. 2 * x * exp (x\<^sup>2)"])
   2.841 -      apply (auto intro!: derivative_eq_intros eventually_elim1[OF eventually_gt_at_top, of 1])
   2.842 -      apply (subst inverse_eq_divide[symmetric])
   2.843 -      apply (rule  tendsto_inverse_0_at_top)
   2.844 -      apply (subst mult_assoc)
   2.845 -      by (auto intro!: filterlim_tendsto_pos_mult_at_top[of "\<lambda>x. 2" 2] filterlim_at_top_mult_at_top[OF filterlim_ident])
   2.846 -    
   2.847 -    then show ?thesis by simp
   2.848 -  qed
   2.849 -
   2.850 -  have "((\<lambda>x. (\<integral>y. (exp (- y\<^sup>2) * indicator {0..} y) * indicator {.. x} y \<partial>lborel) :: real) ---> \<integral>y. exp (- y\<^sup>2) * indicator {0..} y \<partial>lborel) at_top"
   2.851 -    by (intro tendsto_integral_at_top integrable_mult_indicator) auto
   2.852 -  also have "(\<lambda>x. (\<integral>y. (exp (- y\<^sup>2) * indicator {0..} y) * indicator {.. x} y \<partial>lborel) :: real) = 
   2.853 -    (\<lambda>x. (\<integral>y. exp (- y\<^sup>2) * indicator {0..x} y \<partial>lborel) :: real)"
   2.854 -    by (auto simp: fun_eq_iff split: split_indicator intro!: integral_cong)
   2.855 -  finally have *: "((\<lambda>x. (\<integral>y. exp (- y\<^sup>2) * indicator {0..x} y \<partial>lborel) :: real) ---> \<integral>y. exp (- y\<^sup>2) * indicator {0..} y \<partial>lborel) at_top"
   2.856 -    .
   2.857 -
   2.858 -  have tends: "((\<lambda>x. (\<integral> xa. exp (- xa\<^sup>2) * indicator {0..x} xa \<partial>lborel) / 2) ---> (sqrt pi / 2) / 2) at_top"
   2.859 -    apply (rule tendsto_divide)
   2.860 -    apply (subst I[symmetric])
   2.861 -    apply (auto intro: *)
   2.862 -    done
   2.863 -
   2.864 -  have [intro]: "(?IFunc ---> sqrt pi / 4) at_top"
   2.865 -    apply (simp add: byparts)
   2.866 -    apply (subst filterlim_cong[where g = ?f])
   2.867 -    apply (auto simp: eventually_ge_at_top linorder_not_less)
   2.868 -  proof -
   2.869 -    have "((\<lambda>x. (\<integral> xa. exp (- xa\<^sup>2) * indicator {0..x} xa / 2 \<partial>lborel) - x * exp (- x\<^sup>2) / 2::real) --->
   2.870 -        (0 + sqrt pi / 4 - 0)) at_top"
   2.871 -      apply (intro tendsto_diff)
   2.872 -      apply auto
   2.873 -      apply (subst divide_real_def)
   2.874 -      using tends
   2.875 -      by (auto intro!: integrable_mult_indicator)
   2.876 -    then show "((\<lambda>x. (\<integral> xa. exp (- xa\<^sup>2) * indicator {0..x} xa \<partial>lborel) / 2 - x * exp (- x\<^sup>2) / 2) ---> sqrt pi / 4) at_top" by  simp
   2.877 -  qed
   2.878 -    
   2.879 -  have [intro]:"\<And>y. integrable lborel (\<lambda>x. x\<^sup>2 * exp (- x\<^sup>2) * indicator {0..} x * indicator {..y} x::real)"
   2.880 -    apply (subst integrable_cong[where g = "\<lambda>x. x\<^sup>2 * exp (- x\<^sup>2) * indicator {0..y} x" for y])
   2.881 -    by (auto intro!: borel_integrable_atLeastAtMost split:split_indicator)
   2.882 -    
   2.883 -  have **[intro]: "integrable lborel (\<lambda>x. x\<^sup>2 * exp (- x\<^sup>2) * indicator {0..} x::real)"
   2.884 -    by (rule integrable_monotone_convergence_at_top) auto
   2.885 -
   2.886 -  have "(\<integral>x. x\<^sup>2 * exp (- x\<^sup>2) * indicator {0..} x \<partial>lborel) = sqrt pi / 4"
   2.887 -    by (rule integral_monotone_convergence_at_top) auto
   2.888 -  
   2.889 -  then have "(\<integral>\<^sup>+x. ereal (x\<^sup>2 * exp (- x\<^sup>2)* indicator {0..} x) \<partial>lborel) = sqrt pi / 4"
   2.890 -    by (subst nn_integral_eq_integral) auto
   2.891 -
   2.892 -  then have ***: "(\<integral>\<^sup>+ x. ereal (x\<^sup>2 * exp (- x\<^sup>2)) \<partial>lborel) = sqrt pi / 2"
   2.893 -    by (subst nn_integral_even_function)
   2.894 -       (auto simp: real_sqrt_mult real_sqrt_divide ereal_mult_indicator)
   2.895 -  
   2.896 -  then show "(\<integral> x. x\<^sup>2 * exp (- x\<^sup>2) \<partial>lborel) = sqrt pi / 2"
   2.897 -    by (subst integral_eq_nn_integral) auto
   2.898 -
   2.899 -  show "integrable lborel (\<lambda>x. x\<^sup>2 * exp (- x\<^sup>2)::real)"
   2.900 -    by (intro integrableI_nn_integral_finite[OF _ _ ***]) auto
   2.901 -qed
   2.902 -
   2.903 -lemma integral_xsquare_times_standard_normal[intro]: "(\<integral> x. std_normal_density x * x\<^sup>2 \<partial>lborel) = 1"
   2.904 -  and integrable_xsquare_times_standard_normal[intro]: "integrable lborel (\<lambda>x. std_normal_density x * x\<^sup>2)"
   2.905 -proof -
   2.906 -  have [intro]:"integrable lborel (\<lambda>x. exp (- x\<^sup>2) * (2 * x\<^sup>2) / (sqrt 2 * sqrt pi))"
   2.907 -    apply (subst integrable_cong[where g ="(\<lambda>x. (2 * inverse (sqrt 2 * sqrt pi)) * (exp (- x\<^sup>2) * x\<^sup>2))"])
   2.908 -    by (auto intro!: integrable_xsquare_exp_xsquare simp: field_simps)
   2.909 -
   2.910 -  have "(\<integral> x. std_normal_density x * x\<^sup>2 \<partial>lborel) = (2 / sqrt pi) * \<integral> x. x\<^sup>2 * exp (- x\<^sup>2) \<partial>lborel"
   2.911 -    apply (subst integral_mult_right[symmetric])
   2.912 -    apply (rule integrable_xsquare_exp_xsquare)
   2.913 -    unfolding std_normal_density_def
   2.914 -    apply (subst lborel_integral_real_affine[where c = "sqrt 2" and t=0], simp_all)
   2.915 -    unfolding integral_mult_right_zero[symmetric] integral_divide_zero[symmetric]
   2.916 -    apply (intro integral_cong)
   2.917 -    by (auto simp: power_mult_distrib real_sqrt_mult)
   2.918 -  also have "... = 1"
   2.919 -    by (subst integral_xsquare_exp_xsquare, auto)
   2.920 -  finally show "(\<integral> x. std_normal_density x * x\<^sup>2 \<partial>lborel) = 1" .
   2.921 -
   2.922 -  show "integrable lborel (\<lambda>x. std_normal_density x * x\<^sup>2)"
   2.923 -    unfolding std_normal_density_def
   2.924 -    apply (subst lborel_integrable_real_affine_iff[where c = "sqrt 2" and t=0, symmetric])
   2.925 -    by (auto simp: power_mult_distrib real_sqrt_mult)
   2.926 -qed
   2.927 -
   2.928 -lemma 
   2.929 -  assumes [arith]: "0 < \<sigma>"
   2.930 -  shows integral_xsquare_times_normal: "(\<integral> x. normal_density \<mu> \<sigma> x * (x - \<mu>)\<^sup>2 \<partial>lborel) = \<sigma>\<^sup>2"
   2.931 -  and integrable_xsquare_times_normal: "integrable lborel (\<lambda>x. normal_density \<mu> \<sigma> x * (x - \<mu>)\<^sup>2)"
   2.932 -proof -
   2.933 -  have "(\<integral> x. normal_density \<mu> \<sigma> x * (x - \<mu>)\<^sup>2 \<partial>lborel) = \<sigma> * \<sigma> * \<integral> x. std_normal_density x * x\<^sup>2 \<partial>lborel"
   2.934 -    unfolding normal_density_def
   2.935 -    apply (subst lborel_integral_real_affine[ where c = \<sigma> and t = \<mu>])
   2.936 -    apply (auto simp: power_mult_distrib)
   2.937 -    unfolding integral_mult_right_zero[symmetric] integral_divide_zero[symmetric]
   2.938 -    apply (intro integral_cong)
   2.939 -    apply auto
   2.940 -    unfolding normal_density_def
   2.941 -    by (auto simp: real_sqrt_mult field_simps power2_eq_square[symmetric])
   2.942 -    
   2.943 -  also have "... = \<sigma>\<^sup>2" 
   2.944 -    by(auto simp: power2_eq_square[symmetric])
   2.945 -
   2.946 -  finally show "(\<integral> x. normal_density \<mu> \<sigma> x * (x - \<mu>)\<^sup>2 \<partial>lborel) = \<sigma>\<^sup>2" .
   2.947 - 
   2.948 -  show "integrable lborel (\<lambda>x. normal_density \<mu> \<sigma> x * (x - \<mu>)\<^sup>2)"
   2.949 -    unfolding normal_density_def
   2.950 -    apply (subst lborel_integrable_real_affine_iff[ where c = \<sigma> and t = \<mu>, symmetric])
   2.951 -    apply auto
   2.952 -    apply (auto simp: power_mult_distrib)
   2.953 -    apply (subst integrable_cong[where g ="(\<lambda>x. \<sigma> * (std_normal_density x * x\<^sup>2))"], auto)
   2.954 -    unfolding std_normal_density_def
   2.955 -    by (auto simp: field_simps real_sqrt_mult power2_eq_square[symmetric])
   2.956 -qed
   2.957 -  
   2.958 -lemma (in prob_space) standard_normal_distributed_variance:
   2.959 -  fixes a b :: real
   2.960 -  assumes D: "distributed M lborel X std_normal_density"
   2.961 -  shows "variance X = 1"
   2.962 -  apply (subst distributed_integral[OF D, of "(\<lambda>x. (x - expectation X)\<^sup>2)", symmetric])
   2.963 -  by (auto simp: standard_normal_distributed_expectation[OF D])
   2.964 +  using integral_normal_moment_nz_1[OF \<sigma>, of \<mu>] distributed_integral[OF D, of "\<lambda>x. x", symmetric]
   2.965 +  by (auto simp: field_simps)
   2.966  
   2.967  lemma (in prob_space) normal_distributed_variance:
   2.968    fixes a b :: real
   2.969 -  assumes [simp, arith]: " 0 < \<sigma>"
   2.970 +  assumes [simp, arith]: "0 < \<sigma>"
   2.971    assumes D: "distributed M lborel X (normal_density \<mu> \<sigma>)"
   2.972    shows "variance X = \<sigma>\<^sup>2"
   2.973 -proof-  
   2.974 -  have *[intro]: "distributed M lborel (\<lambda>x. (X x - \<mu>) / \<sigma>) (\<lambda>x. ereal (std_normal_density x))"
   2.975 -    by (subst normal_standard_normal_convert[OF assms(1) , symmetric]) fact
   2.976 +  using integral_normal_moment_even[of \<sigma> \<mu> 1]
   2.977 +  by (subst distributed_integral[OF D, symmetric])
   2.978 +     (simp_all add: eval_nat_numeral normal_distributed_expectation[OF assms])
   2.979  
   2.980 -  have "variance X = variance  (\<lambda>x. \<mu> + \<sigma> * ((X x - \<mu>) / \<sigma>) )" by simp
   2.981 -  also have "... = \<sigma>\<^sup>2 * 1"
   2.982 -    apply (subst variance_affine)
   2.983 -    apply (auto intro!: standard_normal_distributed_variance prob_space_normal_density
   2.984 -      simp: distributed_integrable[OF *,of "\<lambda>x. x", symmetric]
   2.985 -      distributed_integrable[OF *,of "\<lambda>x. x\<^sup>2", symmetric] variance_affine
   2.986 -      simp del: integral_divide_zero)
   2.987 -    done
   2.988 -
   2.989 -  finally show ?thesis by simp
   2.990 -qed
   2.991 +lemma (in prob_space) standard_normal_distributed_variance:
   2.992 +  "distributed M lborel X std_normal_density \<Longrightarrow> variance X = 1"
   2.993 +  using normal_distributed_variance[of 1 X 0] by simp
   2.994  
   2.995  lemma (in information_space) entropy_normal_density:
   2.996    assumes [arith]: "0 < \<sigma>"
   2.997 @@ -1446,7 +1645,7 @@
   2.998    also have "\<dots> = (\<integral>x. normal_density \<mu> \<sigma> x * (ln (2 * pi * \<sigma>\<^sup>2)) + (normal_density \<mu> \<sigma> x * (x - \<mu>)\<^sup>2) / \<sigma>\<^sup>2 \<partial>lborel) / (2 * ln b)"
   2.999      by (simp del: minus_add_distrib)
  2.1000    also have "\<dots> = (ln (2 * pi * \<sigma>\<^sup>2) + 1) / (2 * ln b)"
  2.1001 -    by (simp add: integrable_xsquare_times_normal integrable_normal integral_xsquare_times_normal)
  2.1002 +    using integral_normal_moment_even[of \<sigma> \<mu> 1] by (simp add: integrable_normal_moment fact_numeral)
  2.1003    also have "\<dots> = log b (2 * pi * exp 1 * \<sigma>\<^sup>2) / 2"
  2.1004      by (simp add: log_def field_simps ln_mult)
  2.1005    finally show ?thesis .