diff --git a/Mathlib/Analysis/Calculus/ContDiff/FaaDiBruno.lean b/Mathlib/Analysis/Calculus/ContDiff/FaaDiBruno.lean
index bdf8e0685e1056..86447ee09e9046 100644
--- a/Mathlib/Analysis/Calculus/ContDiff/FaaDiBruno.lean
+++ b/Mathlib/Analysis/Calculus/ContDiff/FaaDiBruno.lean
@@ -835,6 +835,26 @@ theorem norm_compAlongOrderedFinpartitionL_le :
     ‖c.compAlongOrderedFinpartitionL 𝕜 E F G‖ ≤ 1 :=
   MultilinearMap.mkContinuousLinear_norm_le _ zero_le_one _
 
+theorem norm_compAlongOrderedFinpartitionL_apply_le (f : F [×c.length]→L[𝕜] G) :
+    ‖c.compAlongOrderedFinpartitionL 𝕜 E F G f‖ ≤ ‖f‖ :=
+  (ContinuousLinearMap.le_of_opNorm_le _ c.norm_compAlongOrderedFinpartitionL_le f).trans_eq
+    (one_mul _)
+
+theorem norm_compAlongOrderedFinpartition_sub_compAlongOrderedFinpartition_le
+    (f₁ f₂ : F [×c.length]→L[𝕜] G) (g₁ g₂ : ∀ i, E [×c.partSize i]→L[𝕜] F) :
+    ‖c.compAlongOrderedFinpartition f₁ g₁ - c.compAlongOrderedFinpartition f₂ g₂‖ ≤
+      ‖f₁‖ * c.length * max ‖g₁‖ ‖g₂‖ ^ (c.length - 1) * ‖g₁ - g₂‖ + ‖f₁ - f₂‖ * ∏ i, ‖g₂ i‖ := calc
+  _ ≤ ‖c.compAlongOrderedFinpartition f₁ g₁ - c.compAlongOrderedFinpartition f₁ g₂‖ +
+      ‖c.compAlongOrderedFinpartition f₁ g₂ - c.compAlongOrderedFinpartition f₂ g₂‖ :=
+    norm_sub_le_norm_sub_add_norm_sub ..
+  _ ≤ ‖f₁‖ * c.length * (max ‖g₁‖ ‖g₂‖) ^ (c.length - 1) * ‖g₁ - g₂‖ + ‖f₁ - f₂‖ * ∏ i, ‖g₂ i‖ := by
+    gcongr ?_ + ?_
+    · refine ((c.compAlongOrderedFinpartitionL 𝕜 E F G f₁).norm_image_sub_le g₁ g₂).trans ?_
+      simp only [Fintype.card_fin]
+      gcongr
+      apply norm_compAlongOrderedFinpartitionL_apply_le
+    · exact c.norm_compAlongOrderedFinpartition_le (f₁ - f₂) g₂
+
 end OrderedFinpartition
 
 /-! ### The Faa di Bruno formula -/
@@ -876,6 +896,79 @@ protected noncomputable def taylorComp
     FormalMultilinearSeries 𝕜 E G :=
   fun n ↦ ∑ c : OrderedFinpartition n, q.compAlongOrderedFinpartition p c
 
+/-- An upper estimate (in terms of `Asymptotics.IsBigO`)
+on the difference between two compositions of Taylor series.
+
+Let `p₁`, `p₂`, `q₁`, `q₂` be four families of formal multilinear series
+depending on a parameter `a`.
+Suppose that the norms of `(p₁ · k)`, `(q₁ · k)`, and `(q₂ · k)` are bounded along a filter `l`
+for all `k ≤ n`.
+Also, suppose that $p₁(a, k) - p₂(a, k) = O(f(a))$, $q₁(a, k) - q₂(a, k) = O(f(a))$
+along `l` for all `k ≤ n`.
+Then the difference between `n`th terms of `(p₁ a).taylorComp (q₁ a)` and `(p₂ a).taylorComp (q₂ a)`
+is `O(f(a))` too.
+
+This lemma can be used, e.g., to show that the composition of two $C^{k+α}$ functions
+is a $C^{k+α}$ function. -/
+theorem taylorComp_sub_taylorComp_isBigO
+    {α H : Type*} [NormedAddCommGroup H] {l : Filter α} {p₁ p₂ : α → FormalMultilinearSeries 𝕜 F G}
+    {q₁ q₂ : α → FormalMultilinearSeries 𝕜 E F} {f : α → H} {n : ℕ}
+    (hp_bdd : ∀ k ≤ n, l.IsBoundedUnder (· ≤ ·) (‖p₁ · k‖))
+    (hpf : ∀ k ≤ n, (fun a ↦ p₁ a k - p₂ a k) =O[l] f)
+    (hq₁_bdd : ∀ k ≤ n, l.IsBoundedUnder (· ≤ ·) (‖q₁ · k‖))
+    (hq₂_bdd : ∀ k ≤ n, l.IsBoundedUnder (· ≤ ·) (‖q₂ · k‖))
+    (hqf : ∀ k ≤ n, (fun a ↦ q₁ a k - q₂ a k) =O[l] f) :
+    (fun a ↦ (p₁ a).taylorComp (q₁ a) n - (p₂ a).taylorComp (q₂ a) n) =O[l] f := by
+  simp only [FormalMultilinearSeries.taylorComp, ← Finset.sum_sub_distrib]
+  refine .sum fun c _ ↦ ?_
+  refine .trans (.of_norm_le fun _ ↦
+    c.norm_compAlongOrderedFinpartition_sub_compAlongOrderedFinpartition_le ..) ?_
+  refine .add ?_ ?_
+  · have H₁ : (p₁ · c.length) =O[l] (1 : α → ℝ) := (hp_bdd _ c.length_le).isBigO_one ℝ
+    have H₂ : ∀ m, (q₁ · (c.partSize m)) =O[l] (1 : α → ℝ) := fun m ↦
+      (hq₁_bdd _ <| c.partSize_le _).isBigO_one ℝ
+    have H₃ : ∀ m, (q₂ · (c.partSize m)) =O[l] (1 : α → ℝ) := fun m ↦
+      (hq₂_bdd _ <| c.partSize_le _).isBigO_one ℝ
+    have H₄ : ∀ m, (fun a ↦ q₁ a (c.partSize m) - q₂ a (c.partSize m)) =O[l] f := fun m ↦
+      hqf _ <| c.partSize_le _
+    rw [← Asymptotics.isBigO_pi] at H₂ H₃ H₄
+    have H₅ := ((H₂.prod_left H₃).norm_left.pow (c.length - 1)).mul H₄.norm_norm
+    simpa [mul_assoc] using H₁.norm_left.mul <| H₅.const_mul_left c.length
+  · have H₁ : (fun a ↦ p₁ a c.length - p₂ a c.length) =O[l] f := hpf _ c.length_le
+    have H₂ : ∀ i, (q₂ · (c.partSize i)) =O[l] (1 : α → ℝ) := fun i ↦
+      (hq₂_bdd _ <| c.partSize_le i).isBigO_one ℝ
+    simpa using H₁.norm_norm.mul <| .finsetProd fun i _ ↦ (H₂ i).norm_left
+
+/-- An upper estimate (in terms of `Asymptotics.IsLittleO`)
+on the difference between two compositions of Taylor series.
+
+Let `p₁`, `p₂`, `q₁`, `q₂` be four families of formal multilinear series
+depending on a parameter `a`.
+Suppose that the norms of `(p₁ · k)`, `(q₁ · k)`, and `(q₂ · k)` are bounded along a filter `l`
+for all `k ≤ n`.
+Also, suppose that $p₁(a, k) - p₂(a, k) = o(f(a))$, $q₁(a, k) - q₂(a, k) = o(f(a))$
+along `l` for all `k ≤ n`.
+Then the difference between `n`th terms of `(p₁ a).taylorComp (q₁ a)` and `(p₂ a).taylorComp (q₂ a)`
+is `o(f(a))` too.
+-/
+theorem taylorComp_sub_taylorComp_isLittleO
+    {α H : Type*} [NormedAddCommGroup H] {l : Filter α} {p₁ p₂ : α → FormalMultilinearSeries 𝕜 F G}
+    {q₁ q₂ : α → FormalMultilinearSeries 𝕜 E F} {f : α → H} {n : ℕ}
+    (hp_bdd : ∀ k ≤ n, l.IsBoundedUnder (· ≤ ·) (‖p₁ · k‖))
+    (hpf : ∀ k ≤ n, (fun a ↦ p₁ a k - p₂ a k) =o[l] f)
+    (hq₁_bdd : ∀ k ≤ n, l.IsBoundedUnder (· ≤ ·) (‖q₁ · k‖))
+    (hq₂_bdd : ∀ k ≤ n, l.IsBoundedUnder (· ≤ ·) (‖q₂ · k‖))
+    (hqf : ∀ k ≤ n, (fun a ↦ q₁ a k - q₂ a k) =o[l] f) :
+    (fun a ↦ (p₁ a).taylorComp (q₁ a) n - (p₂ a).taylorComp (q₂ a) n) =o[l] f := calc
+  _ =O[l] fun a ↦ (fun k : Fin (n + 1) ↦ p₁ a k - p₂ a k,
+                    fun k : Fin (n + 1) ↦ q₁ a k - q₂ a k) := by
+    refine taylorComp_sub_taylorComp_isBigO hp_bdd ?_ hq₁_bdd hq₂_bdd ?_
+    all_goals simp only [← Nat.lt_succ_iff, Nat.forall_lt_iff_fin, ← Asymptotics.isBigO_pi]
+    exacts [Asymptotics.isBigO_fst_prod, Asymptotics.isBigO_snd_prod]
+  _ =o[l] f :=
+    .prod_left (Asymptotics.isLittleO_pi.2 fun k ↦ hpf k (by grind))
+      (Asymptotics.isLittleO_pi.2 fun k ↦ hqf k (by grind))
+
 end FormalMultilinearSeries
 
 theorem analyticOn_taylorComp
diff --git a/Mathlib/Data/Fin/Basic.lean b/Mathlib/Data/Fin/Basic.lean
index 4454e57f17fd6f..59eb274f6b2152 100644
--- a/Mathlib/Data/Fin/Basic.lean
+++ b/Mathlib/Data/Fin/Basic.lean
@@ -33,6 +33,14 @@ open Fin Nat Function
 
 attribute [simp] Fin.succ_ne_zero Fin.castSucc_lt_last
 
+theorem Nat.forall_lt_iff_fin {n : ℕ} {p : ∀ k, k < n → Prop} :
+    (∀ k hk, p k hk) ↔ ∀ k : Fin n, p k k.is_lt :=
+  .symm <| Fin.forall_iff
+
+theorem Nat.exists_lt_iff_fin {n : ℕ} {p : ∀ k, k < n → Prop} :
+    (∃ k hk, p k hk) ↔ ∃ k : Fin n, p k k.is_lt :=
+  .symm <| Fin.exists_iff
+
 /-- Elimination principle for the empty set `Fin 0`, dependent version. -/
 def finZeroElim {α : Fin 0 → Sort*} (x : Fin 0) : α x :=
   x.elim0