We define a Radon transform \(R=R_ l\) from functions on \(\mathbb{P}^ n \mathbb{H}\) to functions on \(Gr(l,n; \mathbb{H})\), the quaternionic Grassmann manifold of all the projective \(l\)-planes in \(\mathbb{P}^ n \mathbb{H}\), by averaging functions on \(\mathbb{P}^ n \mathbb{H}\) over projective \(l\)-planes. Under the assumption \(1 \leq l \leq n - 2\), we show that the range of \(R\) is characterized as a kernel of a certain fourth order invariant differential operator on \(Gr(l,n; \mathbb{H})\). Moreover we investigate the relation between the range-characterizing operator and its radial part.