Consider a function $F(x,y,z): \mathbb{R}^3 \mapsto \mathbb{R}^+$ (i.e., $F(x,y,z) > 0 ~~\forall~ x,y,z$) and consider a set of points in the (3D) space, $\{p_1, p_2, \cdots , p_N\}$.
The problem is: find the straight line $L$ (in the 3D space) that minimizes the moment of inertia $M$, defined as: $$M = \sum_{k=1}^N F(p_k) \cdot D^2(p_k, L)$$ where $D(p, L)$ is the perpendicular distance from point $p$ and line $L$ (i.e., the smallest distance between $p$ and a point on $L$)
Notice that for a constant function (or simply no function, just considering a set of points and minimizing the sum of square perpendicular distances), this would boil down to applying PCA to obtain the first principal component. But I'm not sure how I would handle having a "weight" associated to each point.
You're looking for the weighted PCA formulation. Let's derive it here, for the space $\mathbb R^n$ of $n$ dimensions (in your case, $n=3$). So assume you have $N$ points $\{p_1,\dots,p_N\}$ all in $\mathbb R^n$.
Any line $L$ in $\mathbb R^n$ can be represented by a point $a\in\mathbb R^n$ and a unit directional vector $u\in\mathbb R^n$. In other words, points on a line are parametrized as $a+tu$ where $t\in \mathbb R$.
It's easy to verify that the orthogonal projection of point $p_k$ on that line is given by $m_k=a+\langle u, p_k-a\rangle u$ and the squared distance between the line and $p_k$ is $D^2(p_k, L)=\|p_k-m_k\|^2=\|p_k-a\|^2-\langle u, p_k-a\rangle^2$
Thus, the weighted moment of inertia of the cloud of points w.r.t to $L$ is given by $$I(a,u) = \sum_{k=1}^N F(p_k)\,\cdot\,\left[\|p_k-a\|^2-\langle u, p_k-a\rangle^2\right]\tag{1}$$ First, given $u$, let's minimize this with respect to $a$. The gradient of this function in the $a$ variable is $$\begin{split} \nabla_a I(a,u) &= -2\sum_{k=1}^N F(p_k)\,\cdot\,\left[ (p_k-a)-\langle u,p_k-a\rangle u\right]\\ &=\left(-2\sum_{k=1}^N F(p_k)\right)\left [(a_*-a)+2\langle u, a_*-a\rangle u\right] \end{split}$$ where $a_*$ is the weighted center of mass of the cloud of points $\{p_1,\dots,p_N\}$, $$a_* = \frac{\sum_{k=1}^N F(p_k)\cdot p_k}{\sum_{k=1}^N F(p_k)}$$ Thus that gradient is $0$ for any $a$ such that $a=a_*+tu$ with $t\in\mathbb R$. In particular, the gradient vanishes for $a=a_*$. So we have proved that: $$\boxed{\text{The optimal line $L$ goes through the weighted center of gravity }a_* = \frac{\sum_{k=1}^N F(p_k)\cdot p_k}{\sum_{k=1}^N F(p_k)}}$$ We now plug this back into $(1)$ to obtain $I(a_*, u)$, and minimize it with respect to $u$, under the constraint that $\|u\|^2=1$. Using Lagrange multipliers, we know that there exists $\lambda \in\mathbb R$, such that at the minimum $u_*$, $$\nabla_u I(a_*, u_*)=2\lambda u_*$$ This means $$ \sum_{k=1}^N F(p_k)\,\cdot\,\langle u_*, p_k-a_*\rangle (p_k-a_*)=\lambda u_*\tag{2}$$ Using matrix notations, let $P_C$ be the matrix whose columns are the centered vectors $p_k-a_*$, and let $W$ be the diagonal matrix whose diagonal elements are $F(p_k)$. Then $(2)$ can be rewritten as the following eigenvector problem: $$\boxed{P_CWP_C^T u_*=\lambda u_*}$$ which is usually referred to weighted PCA. The optimal direction for line $L$ is the eigenvector with maximal eigenvalue of this weighted covariance matrix. Why the largest eigenvalue? Because if $u_*$ is a unit eigenvector, then $$\begin{split} I(a_*,u_*) &= \sum_{k=1}^N F(p_k)\|p_k-a_*\|^2 - u_*^TP_CWP_C^T u_*\\ &=\sum_{k=1}^N F(p_k)\|p_k-a_*\|^2 - \lambda \|u_*\|^2\\ &=\underbrace{\sum_{k=1}^N F(p_k)\|p_k-a_*\|^2}_{\text{constant w.r.t. } u_*} -\lambda \end{split}$$ and thus maximizing the eigenvalue is equivalent to minimizing the moment of inertia.
This makes intuitive sense if you assume that the weights are positive integers. In that case, you assume that each point $p_k$ is repeated (duplicated) $F(p_k)$ times in the point cloud, leading to the same matrix formulation as above.