Denote by P-n the space of real algebraic polynomials of degree at most n - 1 and consider a multi-index n := (n(1), ... , n(d)) is an element of N-d, d >= 1, of length vertical bar n vertical bar := n(1) + ... + n(d). Then given the nonnegative weight functions w(j) is an element of L-infinity [a, b], 1 <= j <= d, the polynomial Q is an element of P vertical bar n vertical bar+1 \ {0} is called a multiple orthogonal polynomial relative to n and the weights w(j), 1 <= j <= d, if integral([a,b]) wj(x)x(k)Q(x)d mu = 0, 0 <= k <= n(j) - 1, 1 <= j <= d. The above orthogonality relations are equivalent to the conditions for the L-2 multiple best approximation parallel to Q parallel to(L2(wj)) <= parallel to Q - g parallel to(L2(wj)), for all g is an element of P-nj, 1 <= j <= d. The existence of multiple L-2 orthogonal polynomials easily follows from the solvability of the above linear system. The analogous question for the multiple best L-p approximation, i.e., the existence of an extremal polynomial Q(p) is an element of P vertical bar n vertical bar+1 \ (0) satisfying parallel to Q(p)parallel to(Lp(wj)) <= parallel to Q(p) - g parallel to(Lp(wj)), for all g is an element of P-nj, 1 <= j <= d, poses a more difficult nonlinear problem when 1 <= p <= infinity, p not equal 2. In this paper we shall address this question and verify the existence and uniqueness of multiple L-p orthogonal polynomials under proper conditions. (C) 2013 Elsevier Inc. All rights reserved.