# Theory Continuous_Time_Markov_Chain

```(* Author: Johannes Hölzl <hoelzl@in.tum.de> *)

section ‹Continuous-time Markov chains›

theory Continuous_Time_Markov_Chain
imports Discrete_Time_Markov_Process Discrete_Time_Markov_Chain
begin

subsection ‹Trace Operations: relate @{typ "('a × real) stream"} and @{typ "real ⇒ 'a"}›

partial_function (tailrec) trace_at :: "'a ⇒ (real × 'a) stream ⇒ real ⇒ 'a"
where
"trace_at s ω j = (case ω of (t', s')##ω ⇒ if t' ≤ j then trace_at s' ω j else s)"

lemma trace_at_simp[simp]: "trace_at s ((t', s')##ω) j = (if t' ≤ j then trace_at s' ω j else s)"
by (subst trace_at.simps) simp

lemma trace_at_eq:
"trace_at s ω j = (case sfirst (λx. j < fst (shd x)) ω of ∞ ⇒ undefined | enat i ⇒ (s ## smap snd ω) !! i)"
proof (split enat.split; safe)
assume "sfirst (λx. j < fst (shd x)) ω = ∞"
with sfirst_finite[of "λx. j < fst (shd x)" ω]
have "alw (λx. fst (shd x) ≤ j) ω"
then show "trace_at s ω j = undefined"
by (induction arbitrary: s ω rule: trace_at.fixp_induct) (auto split: stream.split)
next
show "sfirst (λx. j < fst (shd x)) ω = enat n ⟹ trace_at s ω j = (s ## smap snd ω) !! n" for n
proof (induction n arbitrary: s ω)
case 0 then show ?case
by (subst trace_at.simps) (auto simp add: enat_0 sfirst_eq_0 split: stream.split)
next
case (Suc n) show ?case
using sfirst.simps[of "λx. j < fst (shd x)" ω] Suc.prems Suc.IH[of "stl ω" "snd (shd ω)"]
by (cases ω) (auto simp add: eSuc_enat[symmetric] split: stream.split if_split_asm)
qed
qed

lemma trace_at_shift: "trace_at s (smap (λ(t, s'). (t + t', s')) ω) t = trace_at s ω (t - t')"
by (induction arbitrary: s ω rule: trace_at.fixp_induct) (auto split: stream.split)

primcorec merge_at :: "(real × 'a) stream ⇒ real ⇒ (real × 'a) stream ⇒ (real × 'a) stream"
where
"merge_at ω j ω' = (case ω of (t, s) ## ω ⇒ if t ≤ j then (t, s)##merge_at ω j ω' else ω')"

lemma merge_at_simp[simp]: "merge_at (x##ω) j ω' = (if fst x ≤ j then x##merge_at ω j ω' else ω')"
by (cases x) (subst merge_at.code; simp)

subsection ‹Exponential Distribution›

definition exponential :: "real ⇒ real measure"
where
"exponential l = density lborel (exponential_density l)"

lemma space_exponential: "space (exponential l) = UNIV"

lemma sets_exponential[measurable_cong]: "sets (exponential l) = sets borel"

lemma prob_space_exponential: "0 < l ⟹ prob_space (exponential l)"
unfolding exponential_def by (intro prob_space_exponential_density)

lemma AE_exponential: "0 < l ⟹ AE x in exponential l. 0 < x"
unfolding exponential_def using AE_lborel_singleton[of 0] by (auto simp add: AE_density exponential_density_def)

lemma emeasure_exponential_Ioi_cutoff:
assumes "0 < l"
shows "emeasure (exponential l) {x <..} = exp (- (max 0 x) * l)"
proof -
interpret prob_space "exponential l"
unfolding exponential_def using ‹0<l› by (rule prob_space_exponential_density)
have *: "prob {xa ∈ space (exponential l). max 0 x < xa} = exp (- max 0 x * l)"
apply (rule exponential_distributedD_gt[OF _ _ ‹0<l›])
apply (auto simp: exponential_def distributed_def)
apply (subst (6) distr_id[symmetric])
apply (subst (2) distr_cong)
apply simp_all
done
have "emeasure (exponential l) {x <..} = emeasure (exponential l) {max 0 x <..}"
using AE_exponential[OF ‹0<l›] by (intro emeasure_eq_AE) auto
also have "… = exp (- (max 0 x) * l)"
using * unfolding emeasure_eq_measure by (simp add: space_exponential greaterThan_def)
finally show ?thesis .
qed

lemma emeasure_exponential_Ioi:
"0 < l ⟹ 0 ≤ x ⟹ emeasure (exponential l) {x <..} = exp (- x * l)"
using emeasure_exponential_Ioi_cutoff[of l x] by simp

lemma exponential_eq_stretch:
assumes "0 < l"
shows "exponential l = distr (exponential 1) borel (λx. (1/l) * x)"
proof (intro measure_eqI)
fix A assume "A ∈ sets (exponential l)"
then have [measurable]: "A ∈ sets borel"
then have [measurable]: "(λx. x / l) -` A ∈ sets borel"
by (rule measurable_sets_borel[rotated]) simp
have "emeasure (exponential l) A =
(∫⇧+x. ennreal l * (indicator (((*) (1/l) -` A) ∩ {0 ..}) (l * x) * ennreal (exp (- (l * x)))) ∂lborel)"
using ‹0 < l›
by (auto simp: ac_simps emeasure_distr exponential_def emeasure_density exponential_density_def
ennreal_mult zero_le_mult_iff
intro!: nn_integral_cong split: split_indicator)
also have "… = (∫⇧+x. indicator (((*) (1/l) -` A) ∩ {0 ..}) x * ennreal (exp (- x)) ∂lborel)"
using ‹0<l›
apply (subst nn_integral_stretch)
apply (auto simp: nn_integral_cmult)
done
also have "… = emeasure (distr (exponential 1) borel (λx. (1/l) * x)) A"
by (auto simp add: emeasure_distr exponential_def emeasure_density exponential_density_def
intro!: nn_integral_cong split: split_indicator)
finally show "emeasure (exponential l) A = emeasure (distr (exponential 1) borel (λx. (1/l) * x)) A" .

lemma uniform_measure_exponential:
assumes "0 < l" "0 ≤ t"
shows "uniform_measure (exponential l) {t <..} = distr (exponential l) borel ((+) t)" (is "?L = ?R")
proof (rule measure_eqI_lessThan)
fix x
have "0 < emeasure (exponential l) {t<..}"
unfolding emeasure_exponential_Ioi[OF assms] by simp
with assms show "?L {x<..} < ∞"
by (simp add: ennreal_divide_eq_top_iff less_top[symmetric] lessThan_Int_lessThan
emeasure_exponential_Ioi)
have *: "((+) t -` {x<..} ∩ space (exponential l)) = {x - t <..}"
by (auto simp: space_exponential)
show "?L {x<..} = ?R {x<..}"
using assms by (simp add: lessThan_Int_lessThan emeasure_exponential_Ioi divide_ennreal
emeasure_distr * emeasure_exponential_Ioi_cutoff exp_diff[symmetric] field_simps split: split_max)
qed (auto simp: sets_exponential)

lemma emeasure_PiM_exponential_Ioi_finite:
assumes "J ⊆ I" "finite J" "⋀i. i ∈ I ⟹ 0 < R i" "0 ≤ x"
shows "emeasure (Π⇩M i∈I. exponential (R i)) (prod_emb I (λi. exponential (R i)) J (Π⇩E j∈J. {x<..})) = exp (- x * (∑i∈J. R i))"
proof (subst emeasure_PiM_emb)
from assms show "(∏i∈J. emeasure (exponential (R i)) {x<..}) = ennreal (exp (- x * sum R J))"
by (subst prod.cong[OF refl emeasure_exponential_Ioi])
(auto simp add: prod_ennreal exp_sum sum_negf[symmetric] sum_distrib_left)
qed (insert assms, auto intro!: prob_space_exponential)

lemma emeasure_PiM_exponential_Ioi_sequence:
assumes R: "summable R" "⋀i. 0 < R i" "0 ≤ x"
shows "emeasure (Π⇩M i∈UNIV. exponential (R i)) (Π i∈UNIV. {x<..}) = exp (- x * suminf R)"
proof -
let ?R = "λi. exponential (R i)" let ?P = "Π⇩M i∈UNIV. ?R i"
let ?N = "λn::nat. prod_emb UNIV ?R {..<n} (Π⇩E i∈{..<n}. {x<..})"
interpret prob_space ?P
by (intro prob_space_PiM prob_space_exponential R)
have "(Π⇩M i∈UNIV. exponential (R i)) (⋂n. ?N n) = (INF n. (Π⇩M i∈UNIV. exponential (R i)) (?N n))"
by (intro INF_emeasure_decseq[symmetric] decseq_emb_PiE) (auto simp: incseq_def)
also have "… = (INF n. ennreal (exp (- x * (∑i<n. R i))))"
using R by (intro INF_cong emeasure_PiM_exponential_Ioi_finite) auto
also have "… = ennreal (exp (- x * (SUP n. (∑i<n. R i))))"
using R
by (subst continuous_at_Sup_antimono[where f="λr. ennreal (exp (- x * r))"])
(auto intro!: bdd_aboveI2[where M="∑i. R i"] sum_le_suminf summable_mult mult_left_mono
continuous_mult continuous_at_ennreal continuous_within_exp[THEN continuous_within_compose3] continuous_minus
simp: less_imp_le antimono_def image_comp)
also have "… = ennreal (exp (- x * (∑i. R i)))"
using R by (subst suminf_eq_SUP_real) (auto simp: less_imp_le)
also have "(⋂n. ?N n) = (Π i∈UNIV. {x<..})"
by (fastforce simp: prod_emb_def Pi_iff PiE_iff space_exponential)
finally show ?thesis
using R by simp
qed

lemma emeasure_PiM_exponential_Ioi_countable:
assumes R: "J ⊆ I" "countable J" "⋀i. i ∈ I ⟹ 0 < R i" "0 ≤ x" and finite: "integrable (count_space J) R"
shows "emeasure (Π⇩M i∈I. exponential (R i)) (prod_emb I (λi. exponential (R i)) J (Π⇩E j∈J. {x<..})) =
exp (- x * (LINT i|count_space J. R i))"
proof cases
assume "finite J" with assms show ?thesis
by (subst emeasure_PiM_exponential_Ioi_finite)
(auto simp: lebesgue_integral_count_space_finite)
next
assume "infinite J"
let ?R = "λi. exponential (R i)" let ?P = "Π⇩M i∈I. ?R i"
define f where "f = from_nat_into J"
have J_eq: "J = range f" and f: "inj f" "f ∈ UNIV → I"
using from_nat_into_inj_infinite[of J] range_from_nat_into[of J] ‹countable J› ‹infinite J› ‹J ⊆ I›
by (auto simp: inj_on_def f_def simp del: range_from_nat_into)
have Bf: "bij_betw f UNIV J"
unfolding J_eq using inj_on_imp_bij_betw[OF f(1)] .

have summable_R: "summable (λi. R (f i))"
using finite unfolding integrable_bij_count_space[OF Bf, symmetric] integrable_count_space_nat_iff
by (rule summable_norm_cancel)

have "emeasure (Π⇩M i∈UNIV. exponential (R (f i))) (Π i∈UNIV. {x<..}) = exp (- x * (∑i. R (f i)))"
using finite assms unfolding J_eq by (intro emeasure_PiM_exponential_Ioi_sequence[OF summable_R]) auto
also have "(Π⇩M i∈UNIV. exponential (R (f i))) = distr ?P (Π⇩M i∈UNIV. exponential (R (f i))) (λω. λi∈UNIV. ω (f i))"
using R by (intro distr_PiM_reindex[symmetric, OF _ f] prob_space_exponential) auto
also have "… (Π i∈UNIV. {x<..}) = ?P ((λω. λi∈UNIV. ω (f i)) -` (Π i∈UNIV. {x<..}) ∩ space ?P)"
using f(2) by (intro emeasure_distr infprod_in_sets) (auto simp: Pi_iff)
also have "(λω. λi∈UNIV. ω (f i)) -` (Π i∈UNIV. {x<..}) ∩ space ?P = prod_emb I ?R J (Π⇩E j∈J. {x<..})"
by (auto simp: prod_emb_def space_PiM space_exponential Pi_iff J_eq)
also have "(∑i. R (f i)) = (LINT i|count_space J. R i)"
using finite
by (subst integral_count_space_nat[symmetric])
(auto simp: integrable_bij_count_space[OF Bf] integral_bij_count_space[OF Bf])
finally show ?thesis .
qed

lemma AE_PiM_exponential_suminf_infty:
fixes R :: "nat ⇒ real"
assumes R: "⋀n. 0 < R n" and finite: "(∑n. ennreal (1 / R n)) = top"
shows "AE ω in Π⇩M n∈UNIV. exponential (R n). (∑n. ereal (ω n)) = ∞"
proof -
let ?P = "Π⇩M n∈UNIV. exponential (R n)"
interpret prob_space "exponential (R n)" for n
by (intro prob_space_exponential R)
interpret product_prob_space "λn. exponential (R n)" UNIV
proof qed

have AE_pos: "AE ω in ?P. ∀i. 0 < ω i"
unfolding AE_all_countable by (intro AE_PiM_component allI prob_space_exponential R AE_exponential) simp

have indep: "indep_vars (λi. borel) (λi x. x i) UNIV"
using PiM_component
apply (subst P.indep_vars_iff_distr_eq_PiM)
apply (auto simp: restrict_UNIV distr_id2)
apply (subst distr_id2)
apply (intro sets_PiM_cong)
apply (auto simp: sets_exponential cong: distr_cong)
done

have [simp]: "0 ≤ x + x * R i ⟷ 0 ≤ x" for x i
using zero_le_mult_iff[of x "1 + R i"] R[of i] by (simp add: field_simps)

have "(∫⇧+ω. eexp (∑n. - ereal (ω n)) ∂?P) = (∫⇧+ω. (INF n. ∏i<n. eexp (- ereal (ω i))) ∂?P)"
proof (intro nn_integral_cong_AE, use AE_pos in eventually_elim)
fix ω :: "nat ⇒ real" assume ω: "∀i. 0 < ω i"
show "eexp (∑n. - ereal (ω n)) = (⨅n. ∏i<n. eexp (- ereal (ω i)))"
proof (rule LIMSEQ_unique[OF _ LIMSEQ_INF])
show "(λi. ∏i<i. eexp (- ereal (ω i))) ⇢ eexp (∑n. - ereal (ω n))"
using ω by (intro eexp_suminf summable_minus_ereal summable_ereal_pos) (auto intro: less_imp_le)
show "decseq (λn. ∏i<n. eexp (- ereal (ω i)))"
using ω by (auto simp: decseq_def intro!: prod_mono3 intro: less_imp_le)
qed
qed
also have "… = (INF n. (∫⇧+ω. (∏i<n. eexp (- ereal (ω i))) ∂?P))"
proof (intro nn_integral_monotone_convergence_INF_AE')
show "AE ω in ?P. (∏i<Suc n. eexp (- ereal (ω i))) ≤ (∏i<n. eexp (- ereal (ω i)))" for n
using AE_pos
proof eventually_elim
case (elim ω)
show ?case
by (rule prod_mono3) (auto simp: elim le_less)
qed
qed (auto simp: less_top[symmetric])
also have "… = (INF n. (∏i<n. (∫⇧+ω. eexp (- ereal (ω i)) ∂?P)))"
proof (intro INF_cong refl indep_vars_nn_integral)
show "indep_vars (λ_. borel) (λi ω. eexp (- ereal (ω i))) {..<n}" for n
proof (rule indep_vars_compose2[of _ _ _ "λi x. eexp(- ereal x)"])
show "indep_vars (λi. borel) (λi x. x i) {..<n}"
by (rule indep_vars_subset[OF indep]) auto
qed auto
qed auto
also have "… = (INF n. (∏i<n. R i * (∫⇧+x. indicator {0 ..} ((1 + R i) * x) * ennreal (exp (- ((1 + R i) * x))) ∂lborel)))"
by (subst product_nn_integral_component)
(auto simp: field_simps exponential_def nn_integral_density ennreal_mult'[symmetric] ennreal_mult''[symmetric]
exponential_density_def exp_diff exp_minus nn_integral_cmult[symmetric]
intro!: INF_cong prod.cong nn_integral_cong split: split_indicator)
also have "… = (INF n. (∏i<n. ennreal (R i / (1 + R i))))"
proof (intro INF_cong prod.cong refl)
show "R i * (∫⇧+ x. indicator {0..} ((1 + R i) * x) * ennreal (exp (- ((1 + R i) * x))) ∂lborel) =
ennreal (R i / (1 + R i))" for i
using nn_intergal_power_times_exp_Ici[of 0] ‹0 < R i›
by (subst nn_integral_stretch[where c="1 + R i"])
(auto simp: mult.assoc[symmetric] ennreal_mult''[symmetric] less_imp_le mult.commute)
qed
also have "… = (INF n. ennreal (∏i<n. R i / (1 + R i)))"
using R by (intro INF_cong refl prod_ennreal divide_nonneg_nonneg) (auto simp: less_imp_le)
also have "… = (INF n. ennreal (inverse (∏i<n. (1 + R i) / R i)))"
by (subst prod_inversef[symmetric]) simp_all
also have "… = (INF n. inverse (ennreal (∏i<n. (1 + R i) / R i)))"
using R by (subst inverse_ennreal) (auto intro!: prod_pos divide_pos_pos simp: add_pos_pos)
also have "… = inverse (SUP n. ennreal (∏i<n. (1 + R i) / R i))"
by (subst continuous_at_Sup_antimono [where f = inverse])
(auto simp: antimono_def image_comp intro!: continuous_on_imp_continuous_within[OF continuous_on_inverse_ennreal'])
also have "(SUP n. ennreal (∏i<n. (1 + R i) / R i)) = top"
proof (cases "SUP n. ennreal (∏i<n. (1 + R i) / R i)")
case (real r)
have "(λn. ennreal (∏i<n. (1 + R i) / R i)) ⇢ r"
using R unfolding real(2)[symmetric]
by (intro LIMSEQ_SUP monoI ennreal_leI prod_mono2) (auto intro!: divide_nonneg_nonneg add_nonneg_nonneg intro: less_imp_le)
then have "(λn. (∏i<n. (1 + R i) / R i)) ⇢ r"
by (rule tendsto_ennrealD)
(use R real in ‹auto intro!: always_eventually prod_nonneg divide_nonneg_nonneg add_nonneg_nonneg intro: less_imp_le›)
moreover have "(1 + R i) / R i = 1 / R i + 1" for i
using ‹0 < R i› by (auto simp: field_simps)
ultimately have "convergent (λn. ∏i<n. 1 / R i + 1)"
by (auto simp: convergent_def)
then have "summable (λi. 1 / R i)"
using R by (subst summable_iff_convergent_prod) (auto intro: less_imp_le)
moreover have "0 ≤ 1 / R i" for i
using R by (auto simp: less_imp_le)
ultimately show ?thesis
using finite ennreal_suminf_neq_top[of "λi. 1 / R i"] by blast
qed
finally have "(∫⇧+ω. eexp (∑n. - ereal (ω n)) ∂?P) = 0"
by simp
then have "AE ω in ?P. eexp (∑n. - ereal (ω n)) = 0"
by (subst (asm) nn_integral_0_iff_AE) auto
then show ?thesis
using AE_pos
proof eventually_elim
show "(∀i. 0 < ω i) ⟹ eexp (∑n. - ereal (ω n)) = 0 ⟹ (∑n. ereal (ω n)) = ∞" for ω
apply (auto simp del: uminus_ereal.simps simp add: uminus_ereal.simps[symmetric]
intro!: summable_iff_suminf_neq_top intro: less_imp_le)
apply (subst (asm) suminf_minus_ereal)
apply (auto intro!: summable_ereal_pos intro: less_imp_le)
done
qed
qed

subsection ‹Transition Rates›

locale transition_rates =
fixes R :: "'a ⇒ 'a ⇒ real"
assumes R_nonneg[simp]: "⋀x y. 0 ≤ R x y"
assumes R_diagonal_0[simp]: "⋀x. R x x = 0"
assumes finite_weight: "⋀x. (∫⇧+y. R x y ∂count_space UNIV) < ∞"
assumes positive_weight: "⋀x. 0 < (∫⇧+y. R x y ∂count_space UNIV)"
begin

abbreviation S :: "(real × 'a) measure"
where "S ≡ (borel ⨂⇩M count_space UNIV)"

abbreviation T :: "(real × 'a) stream measure"
where "T ≡ stream_space S"

abbreviation I :: "'a ⇒ 'a set"
where "I x ≡ {y. 0 < R x y}"

lemma I_countable: "countable (I x)"
proof -
let ?P = "point_measure UNIV (R x)"
interpret finite_measure ?P
proof
show "emeasure ?P (space ?P) ≠ ∞"
using finite_weight
by (simp add: emeasure_density point_measure_def less_top)
qed
from countable_support emeasure_point_measure_finite2[of "{_}" UNIV "R x"]
show ?thesis
qed

definition escape_rate :: "'a ⇒ real" where
"escape_rate x = ∫y. R x y ∂count_space UNIV"

lemma ennreal_escape_rate: "ennreal (escape_rate x) = (∫⇧+y. R x y ∂count_space UNIV)"
using finite_weight[of x] unfolding escape_rate_def
by (intro nn_integral_eq_integral[symmetric]) (auto simp: integrable_iff_bounded)

lemma escape_rate_pos: "0 < escape_rate x"
using positive_weight unfolding ennreal_escape_rate[symmetric] by simp

lemma nonneg_escape_rate[simp]: "0 ≤ escape_rate x"
using escape_rate_pos[THEN less_imp_le] .

lemma prob_space_exponential_escape_rate: "prob_space (exponential (escape_rate x))"
using escape_rate_pos by (rule prob_space_exponential)

lemma measurable_escape_rate[measurable]: "escape_rate ∈ count_space UNIV →⇩M borel"
by auto

lemma measurable_exponential_escape_rate[measurable]: "(λx. exponential (escape_rate x)) ∈ count_space UNIV →⇩M prob_algebra borel"
by (auto simp: space_prob_algebra sets_exponential prob_space_exponential_escape_rate)

interpretation pmf_as_function .

lift_definition J :: "'a ⇒ 'a pmf" is "λx y.  R x y / escape_rate x"
proof safe
show "0 ≤ R x y / escape_rate x" for x y
by (auto intro!: integral_nonneg_AE divide_nonneg_nonneg R_nonneg simp: escape_rate_def)
show "(∫⇧+y. R x y / escape_rate x ∂count_space UNIV) = 1" for x
using escape_rate_pos[of x]
by (auto simp add: divide_ennreal[symmetric] nn_integral_divide ennreal_escape_rate[symmetric] intro!: ennreal_divide_self)
qed

lemma set_pmf_J: "set_pmf (J x) = I x"
using escape_rate_pos[of x] by (auto simp: set_pmf_iff J.rep_eq less_le)

interpretation exp_esc: pair_prob_space "distr (exponential (escape_rate x)) borel ((+) t)" "J x" for x
proof -
interpret prob_space "distr (exponential (escape_rate x)) borel ((+) t)"
by (intro prob_space.prob_space_distr prob_space_exponential_escape_rate) simp
show "pair_prob_space (distr (exponential (escape_rate x)) borel ((+) t)) (measure_pmf (J x))"
by standard
qed

subsection ‹Continuous-time Kernel›

definition K :: "(real × 'a) ⇒ (real × 'a) measure" where
"K = (λ(t, x). (distr (exponential (escape_rate x)) borel ((+) t)) ⨂⇩M J x)"

interpretation K: discrete_Markov_process "borel ⨂⇩M count_space UNIV" K
proof
show "K ∈ borel ⨂⇩M count_space UNIV →⇩M prob_algebra (borel ⨂⇩M count_space UNIV)"
unfolding K_def
apply measurable
apply (rule measurable_snd[THEN measurable_compose])
apply (auto simp: space_prob_algebra prob_space_measure_pmf)
done
qed

interpretation DTMC: MC_syntax J .

lemma in_space_S[simp]: "x ∈ space S"

lemma in_space_T[simp]: "x ∈ space T"

lemma in_space_lim_stream: "ω ∈ space (K.lim_stream x)"
unfolding K.space_lim_stream space_stream_space[symmetric] by simp

lemma prob_space_K_lim: "prob_space (K.lim_stream x)"
using K.lim_stream[THEN measurable_space] by (simp add: space_prob_algebra)

definition select_first :: "'a ⇒ ('a ⇒ real) ⇒ 'a ⇒ bool"
where "select_first x p y = (y ∈ I x ∧ (∀y'∈I x - {y}. p y < p y'))"

lemma select_firstD1: "select_first x p y ⟹ y ∈ I x"

lemma select_first_unique:
assumes y: "select_first x p y1" " select_first x p y2" shows "y1 = y2"
proof -
have "y1 ≠ y2 ⟹ p y1 < p y2" "y1 ≠ y2 ⟹ p y2 < p y1"
using y by (auto simp: select_first_def)
then show "y1 = y2"
by (rule_tac ccontr) auto
qed

lemma The_select_first[simp]: "select_first x p y ⟹ The (select_first x p) = y"
by (intro the_equality select_first_unique)

lemma select_first_INF:
"select_first x p y ⟹ (INF x∈I x. p x) = p y"
by (intro antisym cINF_greatest cINF_lower bdd_belowI2[where m="p y"])
(auto simp: select_first_def le_less)

lemma measurable_select_first[measurable]:
"(λp. select_first x p y) ∈ (Π⇩M y∈I x. borel) →⇩M count_space UNIV"
using I_countable unfolding select_first_def by (intro measurable_pred_countable pred_intros_conj1') measurable

lemma measurable_THE_select_first[measurable]:
"(λp. The (select_first x p)) ∈ (Π⇩M y∈I x. borel) →⇩M count_space UNIV"
by (rule measurable_THE) (auto intro: select_first_unique I_countable dest: select_firstD1)

lemma sets_S_eq: "sets S = sigma_sets UNIV { {t ..} × A | t A. A ⊆ - I x ∨ (∃s∈I x. A = {s}) }"
proof (subst sets_pair_eq)
let ?CI = "λa::real. {a ..}" let ?Ea = "range ?CI"

show "?Ea ⊆ Pow (space borel)" "sets borel = sigma_sets (space borel) ?Ea"
unfolding borel_Ici by auto
show "?CI`Rats ⊆ ?Ea" "(⋃i∈Rats. ?CI i) = space borel"
using Rats_dense_in_real[of "x - 1" "x" for x] by (auto intro: less_imp_le)

let ?Eb = "Pow (- I x) ∪ (λs. {s}) ` I x"
have "b ∈ sigma_sets UNIV (Pow (- I x) ∪ (λs. {s}) ` I x)" for b
proof -
have "b = (b - I x) ∪ (⋃x∈b ∩ I x. {x})"
by auto
also have "… ∈ sigma UNIV (Pow (- I x) ∪ (λs. {s}) ` I x)"
using I_countable by (intro sets.Un sets.countable_UN') auto
finally show ?thesis
by simp
qed
then show "sets (count_space UNIV) = sigma_sets (space (count_space UNIV)) ?Eb"
by auto
show "countable ({- I x} ∪ (⋃s∈I x. {{s}}))"
using I_countable by auto
show "sets (sigma (space borel × space (count_space UNIV)) {a × b |a b. a ∈ ?Ea ∧ b ∈ ?Eb}) =
sigma_sets UNIV {{t ..} × A |t A. A ⊆ - I x ∨ (∃s∈I x. A = {s})}"
apply simp
apply (intro arg_cong[where f="sigma_sets _"])
apply auto
done
qed (auto intro: countable_rat)

subsection ‹Kernel equals Parallel Choice›

abbreviation PAR :: "'a ⇒ ('a ⇒ real) measure"
where
"PAR x ≡ (Π⇩M y∈I x. exponential (R x y))"

lemma PAR_least:
assumes y: "y ∈ I x"
shows "PAR x {p∈space (PAR x). t ≤ p y ∧ select_first x p y} =
emeasure (exponential (escape_rate x)) {t ..} * ennreal (pmf (J x) y)"
proof -
let ?E = "λy. exponential (R x y)" let ?P' = "Π⇩M y∈I x - {y}. ?E y"
interpret P': prob_space ?P'
by (intro prob_space_PiM prob_space_exponential) simp
have *: "PAR x = (Π⇩M y∈insert y (I x - {y}). ?E y)"
using y by (intro PiM_cong) auto
have "0 < R x y"
using y by simp
have **: "(λ(x, X). X(y := x)) ∈ exponential (R x y) ⨂⇩M Pi⇩M (I x - {y}) (λi. exponential (R x i)) →⇩M PAR x"
using y
apply (subst measurable_cong_sets[OF sets_pair_measure_cong[OF sets_exponential sets_PiM_cong[OF refl sets_exponential]] sets_PiM_cong[OF refl sets_exponential]])
apply measurable
apply (rule measurable_fun_upd[where J="I x - {y}"])
apply auto
done
have "PAR x {p∈space (PAR x). t ≤ p y ∧ (∀y'∈I x-{y}. p y < p y')} =
(∫⇧+ty. indicator {t..} ty * ?P' {p∈space ?P'. ∀y'∈I x-{y}. ty < p y'} ∂?E y)"
unfolding * using ‹y ∈ I x›
apply (subst distr_pair_PiM_eq_PiM[symmetric])
apply (auto intro!: prob_space_exponential simp: emeasure_distr insert_absorb)
apply (subst emeasure_distr[OF **])
subgoal
using I_countable by (auto simp: pred_def[symmetric])
apply (subst P'.emeasure_pair_measure_alt)
subgoal
using I_countable[of x]
apply (intro measurable_sets[OF **])
apply (auto simp: pred_def[symmetric])
done
apply (auto intro!: nn_integral_cong arg_cong2[where f=emeasure] split: split_indicator if_split_asm
simp: space_exponential space_PiM space_pair_measure PiE_iff extensional_def)
done
also have "… = (∫⇧+ty. indicator {t..} ty * ennreal (exp (- ty * (escape_rate x - R x y))) ∂?E y)"
apply (intro nn_integral_cong_AE)
using AE_exponential[OF ‹0 < R x y›]
proof eventually_elim
fix ty :: real assume "0 < ty"
have "escape_rate x =
(∫⇧+y'. R x y' * indicator {y} y' ∂count_space UNIV) + (∫⇧+y'. R x y' * indicator (I x - {y}) y' ∂count_space UNIV)"
unfolding ennreal_escape_rate by (subst nn_integral_add[symmetric]) (auto simp: less_le split: split_indicator intro!: nn_integral_cong)
also have "… = R x y + (∫⇧+y'. R x y' ∂count_space (I x - {y}))"
by (auto simp add: nn_integral_count_space_indicator less_le simp del: nn_integral_indicator_singleton
intro!: arg_cong2[where f="(+)"] nn_integral_cong split: split_indicator)
finally have "(∫⇧+y'. R x y' ∂count_space (I x - {y})) = escape_rate x - R x y ∧ R x y ≤ escape_rate x"
using escape_rate_pos[THEN less_imp_le]
by (cases "(∫⇧+y'. R x y' ∂count_space (I x - {y}))")
(auto simp: add_top ennreal_plus[symmetric] simp del: ennreal_plus)
then have "integrable (count_space (I x - {y})) (R x)" "(LINT y'|count_space (I x - {y}). R x y') = escape_rate x - R x y"
by (auto simp: nn_integral_eq_integrable)
then have "?P' (prod_emb (I x-{y}) ?E (I x-{y}) (Π⇩E j∈(I x-{y}). {ty<..})) = exp (- ty * (escape_rate x - R x y))"
using I_countable ‹0 < ty› by (subst emeasure_PiM_exponential_Ioi_countable) auto
also have "prod_emb (I x-{y}) ?E (I x-{y}) (Π⇩E j∈(I x-{y}). {ty<..}) = {p∈space ?P'. ∀y'∈I x-{y}. ty < p y'}"
by (simp add: set_eq_iff prod_emb_def space_PiM space_exponential ac_simps Pi_iff)
finally show "indicator {t..} ty * ?P' {p∈space ?P'. ∀y'∈I x-{y}. ty < p y'} =
indicator {t..} ty * ennreal (exp (- ty * (escape_rate x - R x y)))"
by simp
qed
also have "… = (∫⇧+ty. ennreal (R x y) * (ennreal (exp (- ty * escape_rate x)) * indicator {max 0 t..} ty) ∂lborel)"
intro!: nn_integral_cong split: split_indicator)
also have "… = (R x y / escape_rate x) * emeasure (exponential (escape_rate x)) {max 0 t..}"
using escape_rate_pos[of x]
by (auto simp: exponential_def exponential_density_def emeasure_density nn_integral_cmult[symmetric] ennreal_mult[symmetric]
split: split_indicator intro!: nn_integral_cong )
also have "… = pmf (J x) y * emeasure (exponential (escape_rate x)) {t..}"
using AE_exponential[OF escape_rate_pos[of x]]
by (intro arg_cong2[where f="(*)"] emeasure_eq_AE) (auto simp: J.rep_eq )
finally show ?thesis
using assms by (simp add: mult_ac select_first_def)
qed

lemma AE_PAR_least: "AE p in PAR x. ∃y∈I x. select_first x p y"
proof -
have D: "disjoint_family_on (λy. {p ∈ space (PAR x). select_first x p y}) (I x)"
by (auto simp: disjoint_family_on_def dest: select_first_unique)
have "PAR x {p∈space (PAR x). ∃y∈I x. select_first x p y} =
PAR x (⋃y∈I x. {p∈space (PAR x). select_first x p y})"
by (auto intro!: arg_cong2[where f=emeasure])
also have "… = (∫⇧+y. PAR x {p∈space (PAR x). select_first x p y} ∂count_space (I x))"
using I_countable by (intro emeasure_UN_countable D) auto
also have "… = (∫⇧+y. PAR x {p∈space (PAR x). 0 ≤ p y ∧ select_first x p y} ∂count_space (I x))"
proof (intro nn_integral_cong emeasure_eq_AE, goal_cases)
case (1 y) with AE_PiM_component[of "I x" "λy. exponential (R x y)" y "(<) 0"] AE_exponential[of "R x y"] show ?case
by (auto simp: prob_space_exponential)
qed (insert I_countable, auto)
also have "… = (∫⇧+y. emeasure (exponential (escape_rate x)) {0 ..} * ennreal (pmf (J x) y) ∂count_space (I x))"
by (auto simp add: PAR_least intro!: nn_integral_cong)
also have "… = (∫⇧+y. emeasure (exponential (escape_rate x)) {0 ..} ∂J x)"
by (auto simp: nn_integral_measure_pmf nn_integral_count_space_indicator ac_simps pmf_eq_0_set_pmf set_pmf_J
simp del: nn_integral_const intro!: nn_integral_cong split: split_indicator)
also have "… = 1"
using AE_exponential[of "escape_rate x"]
by (auto intro!: prob_space.emeasure_eq_1_AE prob_space_exponential simp: escape_rate_pos less_imp_le)
finally show ?thesis
using I_countable
by (subst prob_space.AE_iff_emeasure_eq_1 prob_space_PiM prob_space_exponential)
(auto intro!: prob_space_PiM prob_space_exponential simp del: Set.bex_simps(6))
qed

lemma K_alt: "K (t, x) = distr (Π⇩M y∈I x. exponential (R x y)) S (λp. (t + (INF y∈I x. p y), The (select_first x p)))" (is "_ = ?R")
proof (rule measure_eqI_generator_eq_countable)
let ?E = "{ {t ..} × A | (t::real) A. A ⊆ - I x ∨ (∃s∈I x. A = {s}) }"
show "Int_stable ?E"
apply (auto simp: Int_stable_def)
subgoal for t1 A1 t2 A2
by (intro exI[of _ "max t1 t2"] exI[of _ "A1 ∩ A2"]) auto
subgoal for t1 t2 y1 y2
by (intro exI[of _ "max t1 t2"] exI[of _ "{y1} ∩ {y2}"]) auto
done
show "sets (K (t, x)) = sigma_sets UNIV ?E"
unfolding K.sets_K[OF in_space_S] by (subst sets_S_eq) rule
show "sets ?R = sigma_sets UNIV ?E"
using sets_S_eq by simp
show "countable ((λ(t, A). {t ..} × A) ` (ℚ × ({- I x} ∪ (λs. {s}) ` I x)))"
by (intro countable_image countable_SIGMA countable_rat countable_Un I_countable) auto

have *: "(+) t -` {t'..} ∩ space (exponential (escape_rate x)) = {t' - t..}" for t'
by (auto simp: space_exponential)
{ fix X assume "X ∈ ?E"
then consider
t' s where "s ∈ I x" "X = {t' ..} × {s}"
| t' A where "A ⊆ - I x" "X = {t' ..} × A"
by auto
then show "K (t, x) X = ?R X"
proof cases
case 1
have "AE p in PAR x. (t' - t ≤ p s ∧ select_first x p s) =
(t' ≤ t + (⨅x∈I x. p x) ∧ The (select_first x p) = s)"
using AE_PAR_least by eventually_elim (auto dest: select_first_unique simp: select_first_INF)
with 1 I_countable show ?thesis
by (auto simp add: K_def measure_pmf.emeasure_pair_measure_Times emeasure_distr emeasure_pmf_single *
PAR_least[symmetric] intro!: emeasure_eq_AE)
next
case 2
moreover
then have "emeasure (measure_pmf (J x)) A = 0"
by (subst AE_iff_measurable[symmetric, where P="λx. x ∉ A"])
(auto simp: AE_measure_pmf_iff set_pmf_J subset_eq)
moreover
have "PAR x ((λp. (t + ⨅(p ` (I x)), The (select_first x p))) -` ({t'..} × A) ∩ space (PAR x)) = 0"
using ‹A ⊆ - I x› AE_PAR_least[of x] I_countable
by (subst AE_iff_measurable[symmetric, where P="λp. (t + ⨅(p ` (I x)), The (select_first x p)) ∉ {t'..} × A"])
(auto simp del: all_simps(5) simp add: imp_ex imp_conjL subset_eq)
ultimately show ?thesis
using I_countable
by (simp add: K_def measure_pmf.emeasure_pair_measure_Times emeasure_distr *)
qed }

interpret prob_space "K ts" for ts
by (rule K.prob_space_K) simp
show "emeasure (K (t, x)) a ≠ ∞" for a
using emeasure_finite by simp
qed (insert  Rats_dense_in_real[of "x - 1" x for x], auto, blast intro: less_imp_le)

lemma AE_K: "AE y in K x. fst x < fst y ∧ snd y ∈ J (snd x)"
unfolding K_def split_beta
apply (subst exp_esc.AE_pair_iff[symmetric])
apply measurable
apply (simp_all add: AE_distr_iff AE_measure_pmf_iff exponential_def AE_density exponential_density_def cong del: AE_cong)
using AE_lborel_singleton[of 0]
apply eventually_elim
apply simp
done

lemma AE_lim_stream:
"AE ω in K.lim_stream x. ∀i. snd ((x ## ω) !! i) ∈ DTMC.acc``{snd x} ∧ snd (ω !! i) ∈ J (snd ((x ## ω) !! i)) ∧ fst ((x ## ω) !! i) < fst (ω !! i)"
(is "AE ω in K.lim_stream x. ∀i. ?P ω i")
unfolding AE_all_countable
proof
let ?F = "λi x ω. fst ((x ## ω) !! i)" and ?S = "λi x ω. snd ((x ## ω) !! i)"
fix i show "AE ω in K.lim_stream x. ?P ω i"
proof (induction i arbitrary: x)
case 0 with AE_K[of x] show ?case
by (subst K.AE_lim_stream) (auto simp add: space_pair_measure cong del: AE_cong)
next
case (Suc i)
show ?case
proof (subst K.AE_lim_stream, goal_cases)
case 2 show ?case
using DTMC.countable_reachable
by (intro measurable_compose_countable_restrict[where f="?S (Suc i) x"])
(simp_all del: Image_singleton_iff)
next
case 3 show ?case
apply (simp del: AE_conj_iff cong del: AE_cong)
using AE_K[of x]
apply eventually_elim
subgoal premises K_prems for y
using Suc
by eventually_elim (insert K_prems, auto intro: converse_rtrancl_into_rtrancl)
done
qed
qed

lemma measurable_merge_at[measurable]: "(λ(ω, ω'). merge_at ω j ω') ∈ (T ⨂⇩M T) →⇩M T"
proof (rule measurable_stream_space2)
define F where "F x n = (case x of (ω::(real × 'a) stream, ω') ⇒ merge_at ω j ω') !! n" for x n
fix n
have "(λx. F x n) ∈ stream_space S ⨂⇩M stream_space S →⇩M S"
proof (induction n)
case 0 then show ?case
by (simp add: F_def split_beta' stream.case_eq_if)
next
case (Suc n)
from Suc[measurable]
have eq: "F x (Suc n) = (case fst x of (t, s) ## ω ⇒ if t ≤ j then F (ω, snd x) n else snd x !! Suc n)" for x
by (auto simp: F_def split: prod.split stream.split)
show ?case
unfolding eq stream.case_eq_if by measurable
qed
then show "(λx. (case x of (ω, ω') ⇒ merge_at ω j ω') !! n) ∈ stream_space S ⨂⇩M stream_space S →⇩M S"
unfolding F_def by auto
qed

lemma measurable_trace_at[measurable]: "(λ(s, ω). trace_at s ω j) ∈ (count_space UNIV ⨂⇩M T) →⇩M count_space UNIV"
unfolding trace_at_eq by measurable

lemma measurable_trace_at': "(λ((s, j), ω). trace_at s ω j) ∈ ((count_space UNIV ⨂⇩M borel) ⨂⇩M T) →⇩M count_space UNIV"
unfolding trace_at_eq split_beta' by measurable

lemma K_time_split:
assumes "t ≤ j" and [measurable]: "f ∈ S →⇩M borel"
shows "(∫⇧+x. f x * indicator {j <..} (fst x) ∂K (t, s)) = (∫⇧+x. f x ∂K (j, s)) * exponential (escape_rate s) {j - t <..}"
proof -
have "(∫⇧+ y. ∫⇧+ x. f (t + x, y) * indicator {j<..} (t + x) ∂exponential (escape_rate s) ∂J s) =
(∫⇧+ y. ∫⇧+ x. f (t + x, y) * indicator {j - t<..} x ∂exponential (escape_rate s) ∂J s)"
by (intro nn_integral_cong) (auto split: split_indicator)
also have "… = (∫⇧+ y. ∫⇧+ x. f (t + x, y) ∂uniform_measure (exponential (escape_rate s)) {j-t <..} ∂J s) *
emeasure (exponential (escape_rate s)) {j - t <..}"
using ‹t ≤ j› escape_rate_pos
by (subst nn_integral_uniform_measure)
(auto simp: nn_integral_divide ennreal_divide_times emeasure_exponential_Ioi)
also have "… = (∫⇧+ y. ∫⇧+ x. f (j + x, y) ∂exponential (escape_rate s) ∂J s) *
emeasure (exponential (escape_rate s)) {j - t <..}"
using ‹t ≤ j› escape_rate_pos by (simp add: uniform_measure_exponential nn_integral_distr)
finally show ?thesis
by (simp add: K_def exp_esc.nn_integral_snd[symmetric] nn_integral_distr)
qed

lemma K_in_space[simp]: "K x ∈ space (prob_algebra S)"
by (rule measurable_space [OF K.K]) simp

lemma L_in_space[simp]: "K.lim_stream x ∈ space (prob_algebra T)"
by (rule measurable_space [OF K.lim_stream]) simp

subsection ‹Markov Chain Property›

lemma lim_time_split:
"t ≤ j ⟹ K.lim_stream (t, s) = do { ω ← K.lim_stream (t, s) ; ω' ← K.lim_stream (j, trace_at s ω j) ; return T (merge_at ω j ω')}"
(is "_ ⟹ _ = ?DO t s")
proof (coinduction arbitrary: t s rule: K.lim_stream_eq_coinduct)
case step let ?L = K.lim_stream

note measurable_compose[OF measurable_prob_algebraD measurable_emeasure_subprob_algebra, measurable (raw)]

define B' where "B' = (λ(t', s). if t' ≤ j then ?DO t' s else ?L (t', s))"
show ?case
proof (intro bexI conjI AE_I2)
show [measurable]: "B' ∈ S →⇩M prob_algebra T"
unfolding B'_def by measurable
show "(∃t s. y = (t, s) ∧ B' y = ?DO t s ∧ t ≤ j) ∨ ?L y = B' y" for y
by (cases y; cases "fst y ≤ j") (auto simp: B'_def)
let ?C = "λx. do { ω ← ?L x; ω' ← ?L (j, trace_at s (x##ω) j); return T (merge_at (x##ω) j ω') }"
have "?DO t s = do { x ← K (t, s); ?C x }"
apply (subst K.lim_stream_eq[OF in_space_S])
apply (subst measurable_cong_sets[OF K.sets_K[OF in_space_S] refl])
apply measurable
apply measurable
apply (subst bind_cong[OF refl bind_cong[OF refl bind_return[OF measurable_prob_algebraD]]])
apply measurable
done
also have "… = K (t, s) ⤜ (λy. B' y ⤜ (λω. return T (y ## ω)))" (is "?DO' = ?R")
proof (rule measure_eqI)
have "sets ?DO' = sets T"
by (intro sets_bind'[OF K_in_space]) measurable
moreover have "sets ?R = sets T"
by (intro sets_bind'[OF K_in_space]) measurable
ultimately show "sets ?DO' = sets ?R"
by simp
fix A assume "A ∈ sets ?DO'"
then have A[measurable]: "A ∈ T"
unfolding ‹sets ?DO' = sets T› .
have "?DO' A = (∫⇧+x. ?C x A ∂K (t, s))"
by (subst emeasure_bind_prob_algebra[OF K_in_space]) measurable
also have "… = (∫⇧+x. ?C x A * indicator {.. j} (fst x) ∂K (t, s)) +
(∫⇧+x. ?C x A * indicator {j <..} (fst x) ∂K (t, s))"
by (subst nn_integral_add[symmetric]) (auto intro!: nn_integral_cong split: split_indicator)
also have "(∫⇧+x. ?C x A * indicator {.. j} (fst x) ∂K (t, s)) =
(∫⇧+y. emeasure (B' y ⤜ (λω. return T (y ## ω))) A * indicator {.. j} (fst y) ∂K (t, s))"
proof (intro nn_integral_cong ennreal_mult_right_cong refl arg_cong2[where f=emeasure])
fix x :: "real × 'a" assume "indicator {..j} (fst x) ≠ (0::ennreal)"
then have "fst x ≤ j"
by (auto split: split_indicator_asm)
then show "?C x = (B' x ⤜ (λω. return T (x ## ω)))"
apply (cases x)
apply measurable
apply measurable
apply (subst bind_return)
apply measurable
done
qed
also have "(∫⇧+x. ?C x A * indicator {j <..} (fst x) ∂K (t, s)) =
(∫⇧+y. emeasure (B' y ⤜ (λω. return T (y ## ω))) A * indicator {j <..} (fst y) ∂K (t, s))"
proof -
have *: "(+) t -` {j<..} = {j - t <..}"
by auto

have "(∫⇧+x. ?C x A * indicator {j <..} (fst x) ∂K (t, s)) =
(∫⇧+x. ?L (j, s) A * indicator {j <..} (fst x) ∂K (t, s))"
by (intro nn_integral_cong ennreal_mult_right_cong refl arg_cong2[where f=emeasure])
(auto simp: K.sets_lim_stream bind_return'' bind_const' prob_space_K_lim prob_space_imp_subprob_space split: split_indicator_asm)
also have "… = ?L (j, s) A * exponential (escape_rate s) {j - t <..}"
by (subst nn_integral_cmult) (simp_all add: K_def exp_esc.nn_integral_snd[symmetric] emeasure_distr space_exponential *)
also have "… = (∫⇧+x. emeasure (?L x ⤜ (λω. return T (x ## ω))) A ∂K (j, s)) * exponential (escape_rate s) {j - t <..}"
by (subst K.lim_stream_eq) (auto simp: emeasure_bind_prob_algebra[OF K_in_space _ A])
also have "… = (∫⇧+y. emeasure (?L y ⤜ (λω. return T (y ## ω))) A * indicator {j <..} (fst y) ∂K (t, s))"
using ‹t ≤ j› by (rule K_time_split[symmetric]) measurable
also have "… = (∫⇧+y. emeasure (B' y ⤜ (λω. return T (y ## ω))) A * indicator {j <..} (fst y) ∂K (t, s))"
by (intro nn_integral_cong ennreal_mult_right_cong refl arg_cong2[where f=emeasure])
(auto simp add: B'_def split: split_indicator_asm)
finally show ?thesis .
qed
also have "(∫⇧+y. emeasure (B' y ⤜ (λω. return T (y ## ω))) A * indicator {.. j} (fst y) ∂K (t, s)) +
(∫⇧+y. emeasure (B' y ⤜ (λω. return T (y ## ω))) A * indicator {j <..} (fst y) ∂K (t, s)) =
(∫⇧+y. emeasure (B' y ⤜ (λω. return T (y ## ω))) A ∂K (t, s))"
by (subst nn_integral_add[symmetric]) (auto intro!: nn_integral_cong split: split_indicator)
also have "… = emeasure (K (t, s) ⤜ (λy. B' y ⤜ (λω. return T (y ## ω)))) A"
by (rule emeasure_bind_prob_algebra[symmetric, OF K_in_space _ A]) auto
finally show "?DO' A = emeasure (K (t, s) ⤜ (λy. B' y ⤜ (λω. return T (y ## ω)))) A" .
qed
finally show "?DO t s = K (t, s) ⤜ (λy. B' y ⤜ (λω. return T (y ## ω)))" .
qed

lemma K_eq: "K (t, s) = distr (exponential (escape_rate s) ⨂⇩M J s) S (λ(t', s). (t + t', s))"
proof -
have "distr (exponential (escape_rate s)) borel ((+) t) ⨂⇩M distr (J s) (J s) (λx. x) =
distr (exponential (escape_rate s) ⨂⇩M J s) (borel ⨂⇩M J s) (λ(x, y). (t + x, y))"
proof (intro pair_measure_distr)
interpret prob_space "distr (measure_pmf (J s)) (measure_pmf (J s)) (λx. x)"
by (intro measure_pmf.prob_space_distr) simp
show "sigma_finite_measure (distr (measure_pmf (J s)) (measure_pmf (J s)) (λx. x))"
by unfold_locales
qed auto
also have "… = distr (exponential (escape_rate s) ⨂⇩M J s) S (λ(x, y). (t + x, y))"
by (intro distr_cong refl sets_pair_measure_cong) simp
finally show ?thesis
qed

lemma K_shift: "K (t + t', s) = distr (K (t, s)) S (λ(t, s). (t + t', s))"
unfolding K_eq by (subst distr_distr) (auto simp: comp_def split_beta' ac_simps)

lemma K_not_empty: "space (K x) ≠ {}"
by (simp add: K_def space_pair_measure split: prod.split)

lemma lim_stream_not_empty: "space (K.lim_stream x) ≠ {}"
by (simp add: K.space_lim_stream space_pair_measure split: prod.split)

lemma lim_shift: ― ‹Generalize to bijective function on @{const K.lim_stream} invariant on @{const K}›
"K.lim_stream (t + t', s) = distr (K.lim_stream (t, s)) T (smap (λ(t, s). (t + t', s)))"
(is "_ = ?D t s")
proof (coinduction arbitrary: t s rule: K.lim_stream_eq_coinduct)
case step then show ?case
proof (intro bexI[of _ "λ(t, s). ?D (t - t') s"] conjI)
show "?D t s = K (t + t', s) ⤜ (λy. (case y of (t, s) ⇒ ?D (t - t') s) ⤜ (λω. return T (y ## ω)))"
apply (subst K.lim_stream_eq[OF in_space_S])
apply (subst K_shift)
apply (measurable; fail)
apply (measurable; fail)
apply (subst bind_distr[OF _ measurable_prob_algebraD K_not_empty])
apply (measurable; fail)
apply (measurable; fail)
apply (intro bind_cong refl)
apply (measurable; fail)
apply (measurable; fail)
apply (subst bind_distr[OF _ measurable_prob_algebraD lim_stream_not_empty])
apply (measurable; fail)
apply (measurable; fail)
done
qed (auto cong: conj_cong intro!: exI[of _ "_ - t'"])
qed simp

lemma lim_0: "K.lim_stream (t, s) = distr (K.lim_stream (0, s)) T (smap (λ(t', s). (t' + t, s)))"
using lim_shift[of 0 t s] by simp

subsection ‹Explosion time›

definition explosion :: "(real × 'a) stream ⇒ ereal"
where "explosion ω = (SUP i. ereal (fst (ω !! i)))"

lemma ball_less_Suc_eq: "(∀i<Suc n. P i) ⟷ (P 0 ∧ (∀i<n. P (Suc i)))"
using less_Suc_eq_0_disj by auto

lemma lim_stream_timediff_eq_exponential_1:
"distr (K.lim_stream ts) (PiM UNIV (λ_. borel))
(λω i. escape_rate (snd ((ts##ω) !! i)) * (fst (ω !! i) - fst ((ts##ω) !! i))) =
PiM UNIV (λ_. exponential 1)"
(is "?D = ?P")
proof (rule measure_eqI_PiM_sequence)
show "sets ?D = sets (PiM UNIV (λ_. borel))" "sets ?P = sets (PiM UNIV (λ_. borel))"
by (auto intro!: sets_PiM_cong simp: sets_exponential)
have [measurable]: "ts ∈ space S"
by auto
{ interpret prob_space ?D
by (intro prob_space.prob_space_distr K.prob_space_lim_stream measurable_abs_UNIV) auto
show "finite_measure ?D"
by unfold_locales }

interpret E: prob_space "exponential 1"
by (rule prob_space_exponential) simp
interpret P: product_prob_space "λ_. exponential 1" UNIV
by unfold_locales

let "distr _ _ (?f ts)" = ?D

fix A :: "nat ⇒ real set" and n :: nat assume A[measurable]: "⋀i. A i ∈ sets borel"
define n' where "n' = Suc n"
have "emeasure ?D (prod_emb UNIV (λ_. borel) {..n} (Pi⇩E {..n} A)) =
emeasure (K.lim_stream ts) {ω∈space (stream_space S). ∀i<n'. ?f ts ω i ∈ A i}"
apply (subst emeasure_distr)
apply (auto intro!: measurable_abs_UNIV arg_cong[where f="emeasure _"])
apply (auto simp: prod_emb_def K.space_lim_stream space_pair_measure n'_def)
done
also have "… = (∏i<n'. emeasure (exponential 1) (A i))"
using A
proof (induction n' arbitrary: A ts)
case 0 then show ?case
using prob_space.emeasure_space_1[OF prob_space_K_lim]
next
case (Suc n A ts)
from Suc.prems[measurable]
have [measurable]: "ts ∈ space S"
by auto

have "emeasure (K.lim_stream ts) {ω ∈ space (stream_space S). ∀i<Suc n. ?f ts ω i ∈ A i} =
(∫⇧+ts'. indicator (A 0) (escape_rate (snd ts) * (fst ts' - fst ts)) *
emeasure (K.lim_stream ts') {ω ∈ space (stream_space S). ∀i<n. ?f ts' ω i ∈ A (Suc i)} ∂K ts)"
apply (subst K.emeasure_lim_stream)
apply simp
apply measurable
apply (auto intro!: nn_integral_cong arg_cong2[where f=emeasure] split: split_indicator
simp: ball_less_Suc_eq)
done
also have "… = (∫⇧+ts'. indicator (A 0) (escape_rate (snd ts) * (fst ts' - fst ts)) ∂K ts) *
(∏i<n. emeasure (exponential 1) (A (Suc i)))"
by (subst Suc.IH) (simp_all add: nn_integral_multc)
also have "(∫⇧+ts'. indicator (A 0) (escape_rate (snd ts) * (fst ts' - fst ts)) ∂K ts) =
(∫⇧+t. indicator (A 0) (escape_rate (snd ts) * t) ∂exponential (escape_rate (snd ts)))"
by (simp add: K_def exp_esc.nn_integral_snd[symmetric] nn_integral_distr split: prod.split)
also have "… = emeasure (exponential 1) (A 0)"
using escape_rate_pos[of "snd ts"]
by (subst exponential_eq_stretch) (simp_all add: nn_integral_distr)
also have "emeasure (exponential 1) (A 0) * (∏i<n. emeasure (exponential 1) (A (Suc i))) =
(∏i<Suc n. emeasure (exponential 1) (A i))"
by (rule prod.lessThan_Suc_shift[symmetric])
finally show ?case .
qed
also have "… = emeasure ?P (prod_emb UNIV (λ_. borel) {..<n'} (Pi⇩E {..<n'} A))"
using P.emeasure_PiM_emb[of "{..<n'}" A] by (simp add: prod_emb_def space_exponential)
finally show "emeasure ?D (prod_emb UNIV (λ_. borel) {..n} (Pi⇩E {..n} A)) =
emeasure ?P (prod_emb UNIV (λ_. borel) {..n} (Pi⇩E {..n} A))"
qed

lemma AE_explosion_infty:
assumes bdd: "bdd_above (range escape_rate)"
shows "AE ω in K.lim_stream x. explosion ω = ∞"
proof -
have "escape_rate undefined ≤ (SUP x. escape_rate x)"
using bdd by (intro cSUP_upper) auto
then have SUP_escape_pos: "0 < (SUP x. escape_rate x)"
using escape_rate_pos[of undefined] by simp
then have SUP_escape_nonneg: "0 ≤ (SUP x. escape_rate x)"
by (rule less_imp_le)

have [measurable]: "x ∈ space S" by auto
have "(∑i. 1::ennreal) = top"
by (rule sums_unique[symmetric]) (auto simp: sums_def of_nat_tendsto_top_ennreal)
then have "AE ω in (PiM UNIV (λ_. exponential 1)). (∑i. ereal (ω i)) = ∞"
by (intro AE_PiM_exponential_suminf_infty) auto
then have "AE ω in K.lim_stream x.
(∑i. ereal (escape_rate (snd ((x##ω) !! i)) * (fst (ω !! i) - fst ((x##ω) !! i)))) = ∞"
apply (subst (asm) lim_stream_timediff_eq_exponential_1[symmetric, of x])
apply (subst (asm) AE_distr_iff)
apply (auto intro!: measurable_abs_UNIV)
done
then show ?thesis
using AE_lim_stream
proof eventually_elim
case (elim ω)
then have le: "fst ((x##ω) !! n) ≤ fst ((x ## ω) !! m)" if "n ≤ m" for n m
by (intro lift_Suc_mono_le[OF _ ‹n ≤ m›, of "λi. fst ((x ## ω) !! i)"]) (auto intro: less_imp_le)
have [simp]: "fst x ≤ fst ((x##ω) !! i)" "fst ((x##ω) !! i) ≤ fst (ω !! i)" for i
using le[of "i" "Suc i"] le[of 0 i] by auto

have "(∑i. ereal (escape_rate (snd ((x ## ω) !! i)) * (fst (ω !! i) - fst ((x ## ω) !! i)))) =
(SUP n. ∑i<n. ereal (escape_rate (snd ((x ## ω) !! i)) * (fst (ω !! i) - fst ((x ## ω) !! i))))"
by (intro suminf_ereal_eq_SUP) (auto intro!: mult_nonneg_nonneg)
also have "… ≤ (SUP n. (SUP x. escape_rate x) * (ereal (fst ((x ## ω) !! n)) - ereal (fst x)))"
proof (intro SUP_least SUP_upper2)
fix n
have "(∑i<n. ereal (escape_rate (snd ((x ## ω) !! i)) * (fst (ω !! i) - fst ((x ## ω) !! i)))) ≤
(∑i<n. ereal ((SUP i. escape_rate i) * (fst (ω !! i) - fst ((x ## ω) !! i))))"
using elim bdd by (intro sum_mono) (auto intro!: cSUP_upper)
also have "… = (SUP i. escape_rate i) * (∑i<n. fst ((x ## ω) !! Suc i) - fst ((x ## ω) !! i))"
using elim bdd by (subst sum_ereal) (auto simp: sum_distrib_left)
also have "… = (SUP i. escape_rate i) * (fst ((x ## ω) !! n) - fst x)"
by (subst sum_lessThan_telescope) simp
finally show "(∑i<n. ereal (escape_rate (snd ((x ## ω) !! i)) * (fst (ω !! i) - fst ((x ## ω) !! i))))
≤ (SUP x. escape_rate x) * (```