[proofplan]
We translate involutivity from a condition on vector fields into a condition on the 1-forms $\theta_i$ using Cartan's formula $d\theta(X, Y) = X\theta(Y) - Y\theta(X) - \theta([X,Y])$, which for $X, Y \in V$ simplifies to $d\theta(X, Y) = -\theta([X,Y])$. Completing $\{\theta_1, \ldots, \theta_m\}$ to a local coframe gives dual vector fields, some of which span $V$, and the bracket condition $[X,Y] \in V$ becomes $d\theta_k(X, Y) = 0$ for $k \le m$ and $X, Y$ spanning $V$. Expanding $d\theta_k$ in the full coframe and extracting the coefficient of the spanning 2-forms of $\Lambda^2 V^*$ shows the vanishing condition is exactly equivalent to $d\theta_k$ lying in the ideal generated by $\theta_1, \ldots, \theta_m$ in the exterior algebra, which is the claimed form.
[/proofplan]
[step:Complete $\{\theta_i\}$ to a local coframe and extract the dual basis of vector fields]
Fix $p \in M$. Since $\{\theta_1, \ldots, \theta_m\}$ are linearly independent 1-forms at $p$, by a standard extension argument there exists an open neighbourhood $W$ of $p$ on which they extend to a local coframe
\begin{align*}
\theta_1, \ldots, \theta_m, \theta_{m+1}, \ldots, \theta_n \;\in\; \Omega^1(W),
\end{align*}
linearly independent at every point of $W$. Let $\vartheta_1, \ldots, \vartheta_n \in \Gamma(TW)$ be the dual frame, characterised by
\begin{align*}
\theta_i(\vartheta_j) = \delta_{ij}, \quad 1 \le i, j \le n.
\end{align*}
Then for $1 \le j \le n$,
\begin{align*}
\vartheta_j \in V \iff \theta_i(\vartheta_j) = 0 \text{ for all } 1 \le i \le m \iff \delta_{ij} = 0 \text{ for all } i \le m \iff j > m.
\end{align*}
So $V|_W$ is spanned locally by $\vartheta_{m+1}, \ldots, \vartheta_n$, and $\operatorname{rank}(V) = n - m$.
[/step]
[step:Expand $d\theta_k$ in the coframe basis]
The family $\{\theta_p \wedge \theta_q : 1 \le p < q \le n\}$ is a local basis of $\Lambda^2 T^*W$. Hence for each $k = 1, \ldots, n$ there exist unique smooth functions $a^k_{pq} \in C^\infty(W)$, $p < q$, such that
\begin{align*}
d\theta_k = \sum_{1 \le p < q \le n} a^k_{pq}\, \theta_p \wedge \theta_q.
\end{align*}
Extending $a^k_{pq}$ antisymmetrically by $a^k_{qp} := -a^k_{pq}$ and $a^k_{pp} := 0$, we may write the sum symmetrically as
\begin{align*}
d\theta_k = \frac{1}{2} \sum_{p, q = 1}^n a^k_{pq}\, \theta_p \wedge \theta_q.
\end{align*}
[/step]
[step:Reformulate involutivity using Cartan's formula]
Recall Cartan's formula for the exterior derivative of a 1-form $\theta \in \Omega^1(W)$ acting on vector fields $X, Y \in \Gamma(TW)$:
\begin{align*}
d\theta(X, Y) = X\big(\theta(Y)\big) - Y\big(\theta(X)\big) - \theta([X, Y]).
\end{align*}
The distribution $V$ is involutive iff $[X, Y] \in V$ for all $X, Y \in \Gamma(V)$. Since $\theta_i(X) = \theta_i(Y) = 0$ for all $i \le m$ when $X, Y \in \Gamma(V)$, the first two terms of Cartan's formula vanish for $\theta = \theta_i$, $i \le m$, and
\begin{align*}
d\theta_i(X, Y) = -\theta_i([X, Y]).
\end{align*}
Therefore
\begin{align*}
V \text{ involutive} &\iff \theta_i([X, Y]) = 0 \text{ for all } i \le m, \; X, Y \in \Gamma(V) \\
&\iff d\theta_i(X, Y) = 0 \text{ for all } i \le m, \; X, Y \in \Gamma(V).
\end{align*}
Since $V|_W$ is spanned by $\vartheta_{m+1}, \ldots, \vartheta_n$ and $d\theta_i$ is $C^\infty(W)$-bilinear and antisymmetric, the last condition is equivalent to
\begin{align*}
d\theta_i(\vartheta_j, \vartheta_k) = 0 \quad \text{for all } 1 \le i \le m \text{ and } m < j < k \le n.
\end{align*}
[guided]
The bridge from "vector-field" involutivity to a condition on the 1-forms $\theta_i$ is Cartan's formula, which expresses the exterior derivative of a 1-form $\theta$ in terms of Lie-theoretic data:
\begin{align*}
d\theta(X, Y) = X\theta(Y) - Y\theta(X) - \theta([X, Y]).
\end{align*}
Why is this the right tool? Involutivity is about brackets, while $d$ is about antisymmetric derivation. Cartan's formula relates the two, and the exchange becomes very clean when $\theta$ annihilates both $X$ and $Y$.
Specifically, take $\theta = \theta_i$ with $i \le m$ and $X, Y \in \Gamma(V)$. By definition of $V$, $\theta_i(X) = \theta_i(Y) = 0$, so $X\theta_i(Y) = X(0) = 0$ and likewise for the middle term. The formula collapses to
\begin{align*}
d\theta_i(X, Y) = -\theta_i([X, Y]).
\end{align*}
The right side vanishes iff $[X, Y]$ is annihilated by $\theta_i$, and "annihilated by every $\theta_i$, $i \le m$" is exactly the membership condition for $V$. Ranging over all $i \le m$ gives $[X, Y] \in V$ iff $d\theta_i(X, Y) = 0$ for all $i \le m$.
To turn this pointwise condition into a finite test, use that $V$ is spanned locally by the dual vector fields $\vartheta_{m+1}, \ldots, \vartheta_n$. Since $d\theta_i$ is $C^\infty(W)$-bilinear on vector fields (it is an element of $\Omega^2$), vanishing on all of $\Gamma(V) \times \Gamma(V)$ is equivalent to vanishing on basis pairs $(\vartheta_j, \vartheta_k)$ with $m < j < k \le n$. Hence:
\begin{align*}
V \text{ involutive} \iff d\theta_i(\vartheta_j, \vartheta_k) = 0 \text{ for all } 1 \le i \le m, \; m < j < k \le n.
\end{align*}
[/guided]
[/step]
[step:Extract coefficients: identify $d\theta_i(\vartheta_j, \vartheta_k)$ with $a^i_{jk}$]
From Step 2, for $1 \le j < k \le n$,
\begin{align*}
d\theta_i(\vartheta_j, \vartheta_k) = \sum_{1 \le p < q \le n} a^i_{pq}\, (\theta_p \wedge \theta_q)(\vartheta_j, \vartheta_k).
\end{align*}
Using $(\theta_p \wedge \theta_q)(\vartheta_j, \vartheta_k) = \theta_p(\vartheta_j) \theta_q(\vartheta_k) - \theta_p(\vartheta_k) \theta_q(\vartheta_j) = \delta_{pj}\delta_{qk} - \delta_{pk}\delta_{qj}$, only the term $(p, q) = (j, k)$ (which has $p < q$ since $j < k$) contributes, giving
\begin{align*}
d\theta_i(\vartheta_j, \vartheta_k) = a^i_{jk}.
\end{align*}
Combining with the previous step,
\begin{align*}
V \text{ involutive} \iff a^i_{jk} = 0 \text{ for all } 1 \le i \le m \text{ and } m < j < k \le n.
\end{align*}
[/step]
[step:Translate the coefficient condition into the ideal condition $d\theta_i \in (\theta_1, \ldots, \theta_m)$]
We prove that
\begin{align*}
a^i_{jk} = 0 \text{ for all } m < j < k \le n \iff d\theta_i = \sum_{j=1}^m \theta_j \wedge \alpha_{ji} \text{ for some } \alpha_{ji} \in \Omega^1(W).
\end{align*}
$(\Leftarrow)$ Suppose $d\theta_i = \sum_{j=1}^m \theta_j \wedge \alpha_{ji}$. Expand each $\alpha_{ji} = \sum_{\ell=1}^n c_{ji\ell}\, \theta_\ell$ for some $c_{ji\ell} \in C^\infty(W)$. Then
\begin{align*}
d\theta_i = \sum_{j=1}^m \sum_{\ell=1}^n c_{ji\ell}\, \theta_j \wedge \theta_\ell.
\end{align*}
Every term has at least one index ($j$) in the range $\{1, \ldots, m\}$. Comparing with the expansion $d\theta_i = \sum_{p < q} a^i_{pq}\, \theta_p \wedge \theta_q$ and noting the $\theta_p \wedge \theta_q$, $p < q$, are linearly independent, we deduce $a^i_{pq} = 0$ whenever $p, q > m$. Specialising to $m < j < k \le n$ gives $a^i_{jk} = 0$.
$(\Rightarrow)$ Conversely, assume $a^i_{jk} = 0$ whenever $m < j < k \le n$. Then in the expansion
\begin{align*}
d\theta_i = \sum_{1 \le p < q \le n} a^i_{pq}\, \theta_p \wedge \theta_q = \sum_{\substack{1 \le p < q \le n \\ p \le m}} a^i_{pq}\, \theta_p \wedge \theta_q,
\end{align*}
every surviving term has $p \le m$. Rewriting and collecting by the "small" index $p$:
\begin{align*}
d\theta_i = \sum_{p=1}^m \theta_p \wedge \Big(\sum_{q = p+1}^n a^i_{pq}\, \theta_q\Big) = \sum_{p=1}^m \theta_p \wedge \alpha_{pi}, \quad \text{where } \alpha_{pi} := \sum_{q = p+1}^n a^i_{pq}\, \theta_q.
\end{align*}
Each $\alpha_{pi} \in \Omega^1(W)$, and the identity $d\theta_i = \sum_{j=1}^m \theta_j \wedge \alpha_{ji}$ holds with the renaming $j = p$.
[guided]
This step is the algebraic heart of the theorem. We want to translate the coefficient vanishing $a^i_{jk} = 0$ for $j, k > m$ into the statement that $d\theta_i$ lies in the ideal of $\Omega^*(W)$ generated by $\theta_1, \ldots, \theta_m$.
$(\Leftarrow)$ Assume the ideal representation $d\theta_i = \sum_{j=1}^m \theta_j \wedge \alpha_{ji}$. Why does this force the coefficients $a^i_{jk}$ for $j, k > m$ to vanish? Expand each $\alpha_{ji}$ in the coframe: $\alpha_{ji} = \sum_{\ell = 1}^n c_{ji\ell} \theta_\ell$. Substituting,
\begin{align*}
d\theta_i = \sum_{j=1}^m \sum_{\ell=1}^n c_{ji\ell}\, \theta_j \wedge \theta_\ell.
\end{align*}
Every wedge product in this expansion involves at least one factor $\theta_j$ with $j \le m$. Equating with the basis expansion $d\theta_i = \sum_{p<q} a^i_{pq} \theta_p \wedge \theta_q$ and using linear independence of the basis $\{\theta_p \wedge \theta_q\}_{p<q}$, we conclude $a^i_{pq} = 0$ whenever $p, q$ are both $> m$.
$(\Rightarrow)$ Conversely, assume $a^i_{jk} = 0$ for $m < j < k \le n$. Then the expansion of $d\theta_i$ has no "pure $V^*$" terms:
\begin{align*}
d\theta_i = \sum_{\substack{p < q \\ p \le m}} a^i_{pq}\, \theta_p \wedge \theta_q.
\end{align*}
In every surviving term, $\theta_p$ with $p \le m$ sits on the left, so we can pull it out:
\begin{align*}
d\theta_i = \sum_{p=1}^m \theta_p \wedge \Big(\sum_{q > p} a^i_{pq}\, \theta_q\Big) =: \sum_{p=1}^m \theta_p \wedge \alpha_{pi}.
\end{align*}
The "witness" 1-forms $\alpha_{pi} := \sum_{q > p} a^i_{pq} \theta_q$ are smooth because the $a^i_{pq}$ are and the $\theta_q$ are. Renaming the index $p$ to $j$ gives the required ideal representation.
Geometrically, the ideal $I(V) := (\theta_1, \ldots, \theta_m) \subseteq \Omega^*(W)$ is the annihilator ideal of $V$ in the exterior algebra: a form lies in $I(V)$ iff its restriction to $V$ (viewed fibrewise) vanishes. What we have shown is that $V$ is involutive iff the ideal $I(V)$ is closed under $d$, i.e. $d(I(V)) \subseteq I(V)$ — this is the content of the Frobenius involutivity criterion in its differential-ideal form.
[/guided]
[/step]
[step:Combine the equivalences to conclude]
Chaining the equivalences from Steps 3, 4, and 5:
\begin{align*}
V \text{ involutive} &\iff d\theta_i(\vartheta_j, \vartheta_k) = 0 \text{ for all } 1 \le i \le m, \; m < j < k \le n \quad (\text{Step 3}) \\
&\iff a^i_{jk} = 0 \text{ for all } 1 \le i \le m, \; m < j < k \le n \quad (\text{Step 4}) \\
&\iff \exists \alpha_{ji} \in \Omega^1(W) \text{ with } d\theta_i = \sum_{j=1}^m \theta_j \wedge \alpha_{ji} \text{ for each } 1 \le i \le m \quad (\text{Step 5}).
\end{align*}
The 1-forms $\alpha_{ji}$ are defined locally on the chosen chart $W$, but involutivity is a local condition — it may be checked on any open cover — so the equivalence descends to a statement on $M$. This completes the proof.
[/step]