Orthogonal polynomials with respect to weighted inner product Recently I have posted a question on SO, but maybe here is better place to ask.
So, I have data and I want to fit polynomial of order $k$ orthogonal with respect to weighted inner product of functions:
$$
\left<f, g\right> = \sum_{j}^{n} \omega_j f(x_j) g(x_j)
$$
Here $w_j$ is a vector of weights (constants assigned to each of $n$ points $x_j$).
When all $w_j = 1$ we have ordinary euclidean inner product $\left<f, g\right> = \sum f(x_j) g(x_j)$ and we can do it by QR decomposition of Vandermonde matrix $V_{ij}$:
$$ V = \begin{bmatrix} 
1 & x_1 & x_1^2 &... & x_1^k \\ 
1 & x_2 & x_2^2 &... & x_2^k\\ 
... & ... & ... & ...\\
1 & x_n & x_n^2 & ... & x_n^k
\end{bmatrix}.$$
QR decomposition:
$$
V = QR, \quad QQ^T = 1.
$$
So the columns of $Q$ are orthogonal to each other with respect to ordinary euclidean inner product $\left<f, g\right> = \sum f(x_j) g(x_j)$.
My question is: How can I use QR decomposition to get the matrix, whose columns are orthogonal  with respect to weighted inner product of functions?
I was thinking about introducing diagonal weights matrix
$$ W = \begin{bmatrix} 
\omega_1 & 0 & ...\\ 
0 & \omega_2 & ...\\ 
... & ... & ...
\end{bmatrix}$$
and applying QR decomposition to matrix $\sqrt{W}V$:
$$
\sqrt{W}V = \widetilde{Q}R, \quad \widetilde{Q}\widetilde{Q}^T = 1
$$
and after that defining Q as 
$$
Q = \left(\sqrt{W}\right)^{-1}\widetilde{Q},
$$
but I haven't succeeded proving that and maybe it's wrong.
Could anyone help me? I know how to generate such polynomial recursively, but I want to know if one could somehow use QR decomposition.
 A: Not sure about the QR decomposition (see edit), but I found a way to do it with an eigendecomposition (derivation below).
Formulating the problem
Let $n \times (k+1)$ matrix $V$ be the Vandermonde matrix, as defined in the original question. We want to factor $V$ into a product of two matrices: 1) An $n \times (k+1)$ matrix $A$, where $A_{ij}$ contains the value of the $j$th orthogonal polynomial evaluated at the $i$th data point, and 2) A $(k+1) \times (k+1)$ matrix $B$:
$$A B = V$$
The columns of $A$ must be orthogonal in the sense defined in the original question. That is:
$$A^T W A = I$$
where $I$ is the identity matrix and $W$ is a diagonal $n \times n$ matrix containing the weights. The solution below actually applies for any symmetric, positive definite $W$; it doesn't have to be diagonal.
Solution
One possible solution is:
$$A = V U \Lambda^{-\frac{1}{2}} \quad
B = \Lambda^{\frac{1}{2}} U^T$$
where $U \Lambda U^T$ is the eigendecomposition of $V^T W V$
The solution is actually only specified up to unitary transformations (i.e. rotations and sign flips), so there are infinitely many. That is, given a solution $(A,B)$ and any matrix $M$ where $M^T M = I$, then $(A M, M^T B)$ is also a solution.
Derivation
Since $AB = V$ we can write:
$$(A B)^T W (A B)
\enspace = \enspace B^T A^T W A B
\enspace = \enspace V^T W V$$
Since $A^T W A = I$ this simplifies to:
$$B^T B = V^T W V$$
Let $U$ be an orthogonal matrix containing the eigenvectors of $V^T W V$ on the columns, and let diagonal matrix $\Lambda$ contain the corresponding eigenvalues. So $U \Lambda U^T = V^T W V$. Then a solution to the above equation is:
$$B = \Lambda^\frac{1}{2} U^T$$
Plugging this back into $A B = V$ gives:
$$A \Lambda^{\frac{1}{2}} U^T = V$$
Since $U$ is orthogonal and $\Lambda$ is diagonal, we can solve for $A$ by right-multiplying both sides by $U \Lambda^{-\frac{1}{2}}$:
$$A = V U \Lambda^{-\frac{1}{2}}$$
Edit: Solution based on QR decomposition
The OP proposed the following solution (I've changed some variable names): Let $\tilde{Q} \tilde{R}$ be the QR decomposition of $W^{\frac{1}{2}} V$. Then:
$$A = W^{-\frac{1}{2}} \tilde{Q}$$
Here's a proof that the columns are orthogonal:
$$A^T W A
\tag{This must equal the identity matrix}$$
$$= \tilde{Q}^T (W^{-\frac{1}{2}})^T W W^{-\frac{1}{2}} \tilde{Q}
\tag{Substitute in expression for $A$}$$
$$= \tilde{Q}^T W^{-\frac{1}{2}} W W^{-\frac{1}{2}} \tilde{Q}
\tag{$W$ is symmetric}$$
$$= \tilde{Q}^T \tilde{Q}
\tag{$W^{-\frac{1}{2}} W W^{-\frac{1}{2}} = I$}$$
$$= I
\tag{By definition of the QR decomposition}$$
We can find the corresponding $B$ by solving $A B = V$, which gives:
$$B = \tilde{Q}^T W^{\frac{1}{2}} V$$
Therefore, this is also a valid solution. The solution above (based on the eigendecomposition; call it $A_{eig}$) and the OP's solution (based on the QR decomposition; call it $A_{QR}$) are related as:
$$A_{eig} = A_{QR} M
\quad B_{eig} = M^T B_{QR}
\quad \text{where} \quad
M = \tilde{Q}^T W^{\frac{1}{2}} V U L^{-\frac{1}{2}}$$
One can show that $M^T M = I$, so this is a unitary transformation (see the note above about infinitely many solutions).
