A function \(f: V^m\to W\) (with \(V\) and \(W\) being vector spaces over \(\mathbb{Q}\) and a positive integer \(m\)) is called \textit{multi-Jensen} if for each of its \(m\) arguments it satisfies the Jensen equation \(F\left(\frac{x+y}{2}\right)=\frac{F(x)+F(y)}{2}\). Multi-Jensen functions were considered in details by the authors in [Aequationes Math. 69, No. 1--2, 41--57 (2005; Zbl 1072.39025)]. The present paper is devoted to the stability problem. First, it is proved that multi-Jensen mappings can be defined equivalently as solutions of one functional equation: \[ f\left(\frac{x+y}{2}\right)=\frac{1}{2^m}\sum_{S\subseteq\mathbf{m}}f \left(x_S+y_{\mathbf{m}\setminus S}\right),\qquad x,y\in V^m \] where \(\mathbf{m}=\{1,\ldots, m\}\), \(x_S=(0,\ldots,0,x_{j_1},0,\dots,0,x_{j_i},0,\ldots,0)\) for \(S=\{j_1,\ldots,j_i\}\subseteq \mathbf{m}\) and \(x=(x_1,\dots,x_m)\). Assuming that the target space \(W\) is a Banach space, authors prove that if \(f\) is an \(\varepsilon\)-multi-Jensen function, i.e., \[ \left\| f\left(\frac{x+y}{2}\right)-\frac{1}{2^m}\sum_{S\subseteq\mathbf{m}}f \left(x_S+y_{\mathbf{m}\setminus S}\right)\right\| \leq \varepsilon,\qquad x,y\in V^m, \] then there exists a unique (up to a constant) multi-Jensen function \(g: V^m\to W\) such that \(\| f(x)-g(x)\| \leq 2m\varepsilon\) for all \(x\in V^m\). The laborious proof is based on the standard, so called direct, method.