Let $L = L_0+V$ be the higher order Schrödinger type operator where $L_0$ is a homogeneous elliptic operator of order $2m$ in divergence form with bounded coefficients and $V$ is a real measurable function as multiplication operator (e.g., including $(−∆) ^m+V (m∈\mathbb{N})$ as special examples). In this paper, assume that $V$ satisfies a strongly subcritical form condition associated with $L_0$, the authors attempt to establish a theory of Hardy space $H^p_L(\mathbb{R}^n) (0 < p ≤ 1)$ associated with the higher order Schrödinger type operator $L$. Specifically, we first define the molecular Hardy space $H^p_L(\mathbb{R}^n)$ by the so-called $(p,q,ε,M)$ molecule associated to $L$ and then establish its characterizations by the area integral defined by the heat semigroup $e^{−tL}$.