Okay, I simplified my solution. My original solution was about 3 times longer than this, so it is so lucky for me to find a shortcut like this :)
Step 1. Let $I$ be the integral in question:
\begin{align*}
I &= \int_{0}^{1} \frac{\log x \log (1 - x) \log^{2} (1 + x)}{x} \, dx.
\end{align*} By the simple algebraic formula $(a + b)^{3} + (a - b)^{3} - 2 a^{3} = 6 a b^{2}$, it follows that
\begin{align*}
I
&= \frac{1}{6} \int_{0}^{1} \frac{\log x}{x} \left\{ \log^{3} (1-x^{2}) + \log^{3} \left( \frac{1-x}{1+x} \right) - 2\log^{3}(1-x) \right\} \, dx \\
&= \frac{1}{6} \int_{0}^{1} \frac{\log x \log^{3} (1-x^{2})}{x} \, dx + \frac{1}{6} \int_{0}^{1} \frac{\log x \log^{3} \left( \frac{1-x}{1+x} \right)}{x} \, dx - \frac{1}{3} \int_{0}^{1} \frac{\log x \log^{3} (1-x)}{x} \, dx \tag{1}
\end{align*} Applying the substitution $x^{2} \mapsto x$, the first integral reduces to
\begin{align*}
\frac{1}{6} \int_{0}^{1} \frac{\log x \log^{3} (1-x^{2})}{x} \, dx
&= \frac{1}{24} \int_{0}^{1} \frac{\log x \log^{3} (1-x)}{x} \, dx
= \frac{1}{24} \frac{\partial^{4}\beta}{\partial z \partial w^{3}}(0^{+}, 1) \\
&= \frac{1}{2}\zeta(5) - \frac{1}{4}\zeta(2)\zeta(3).
\end{align*} So $\text{(1)}$ can be written as
\begin{align*}
I &= \frac{7}{4} \zeta(2)\zeta(3) - \frac{7}{2} \zeta(5) + I_{2}, \tag{2}
\end{align*} where ( I_{2} ) is given by
\begin{align*}
I_{2} = \frac{1}{6} \int_{0}^{1} \frac{\log x \log^{3} \left( \frac{1-x}{1+x} \right)}{x} \, dx
\end{align*}
Step 2. Now we calculate $I_{2}$. Applying integrating by parts, followed by the substitution $x \mapsto \displaystyle \frac{1-y}{1+y}$, we have
\begin{align*}
I_{2} &= \left[ \frac{1}{12} \log^{2} x \log^{3} \left( \frac{1-x}{1+x} \right) \right]_{0}^{1} + \frac{1}{2} \int_{0}^{1} \frac{\log^{2} x \log^{2} \left( \frac{1-x}{1+x} \right)}{1 - x^{2}} \, dx \\
&= \frac{1}{4} \int_{0}^{1} \frac{\log^{2} y \log^{2} \left( \frac{1-y}{1+y} \right)}{y} \, dy \\
&= \int_{0}^{1} \frac{\log^{2} y}{y} \left\{ \frac{1}{2} \log^{2} \left( \frac{1+y}{1-y} \right) \right\}^{2} \, dy \\
&= \sum_{m=0}^{\infty} \sum_{n=0}^{\infty} \frac{1}{(2m+1)(2n+1)} \int_{0}^{1} y^{2m+2n+1} \log^{2} y \, dy \\
&= \sum_{m=0}^{\infty} \sum_{n=0}^{\infty} \frac{2}{(2m+1)(2n+1)(2m+2n+2)^{3}}.
\end{align*} By the following partial fraction decomposition
\begin{align*}
\frac{1}{ab(a+b)^{3}} = \frac{a+b}{a b(a+b)^{4}} = \frac{1}{a(a+b)^{4}} + \frac{1}{b(a+b)^{4}},
\end{align*}
together with the symmetry in the role of $m$ and $n$, it follows that
\begin{align*}
I_{2}
&= 4 \sum_{m=0}^{\infty} \sum_{n=0}^{\infty} \frac{1}{(2n+1)(2m+2n+2)^{4}}
= \frac{1}{4} \sum_{m=0}^{\infty} \sum_{n=0}^{\infty} \frac{1}{(2n+1)(m+n+1)^{4}} \\
&= \frac{1}{4} \sum_{n=0}^{\infty} \sum_{m=n+1}^{\infty} \frac{1}{2n+1} \cdot \frac{1}{m^{4}}
= \frac{1}{4} \sum_{n=1}^{\infty} \sum_{m=n}^{\infty} \frac{1}{2n-1} \cdot \frac{1}{m^{4}} \\
&= \frac{1}{4} \sum_{m=1}^{\infty} \left( \sum_{n=1}^{m} \frac{1}{2n-1} \right) \frac{1}{m^{4}}.
\end{align*} Now by exploiting the following identity
\begin{align*}
&\int_{0}^{1} \frac{1 - x^{k}}{1 - x} \, dx = 1 + \frac{1}{2} + \cdots + \frac{1}{k} \\
\Longrightarrow & \quad
\int_{0}^{1} \frac{1 + x^{k} - 2x^{2k}}{2(1 - x)} \, dx = 1 + \frac{1}{3} + \cdots + \frac{1}{2k-1},
\end{align*} we can rephrase $I_{2}$ as follows:
\begin{align*}
I_{2}
&= \frac{1}{8} \sum_{m=1}^{\infty} \frac{1}{m^{4}} \int_{0}^{1} \frac{1 + x^{m} - 2x^{2m}}{1 - x} \, dx \\
&= \frac{1}{8} \int_{0}^{1} \frac{1}{1-x} \left( \zeta(4) + \mathrm{Li}_{4}(x) - 2 \mathrm{Li}_{4}(x^{2}) \right) \, dx.
\end{align*} Writing ( \displaystyle \frac{1}{1-x} = \frac{\mathrm{Li}{0}(x)}{x} ) and applying the simple differentiation rule $$\displaystyle \frac{d}{dx} \mathrm{Li}_{s} (x^{p}) = \frac{p \mathrm{Li}_{s-1}(x^{p})}{x}$$, by applying integrating parts 3 times,
\begin{align*}
I_{2}
&= \frac{1}{8} \left[ \mathrm{Li}_{1}(x) \left( \zeta(4) + \mathrm{Li}_{4}(x) - 2 \mathrm{Li}_{4}(x^{2}) \right) \right]_{0}^{1} - \frac{1}{8} \int_{0}^{1} \frac{\mathrm{Li}_{1}(x)}{x} \left( \mathrm{Li}_{3}(x) - 4 \mathrm{Li}_{3}(x^{2}) \right) \, dx \\
&= -\frac{1}{8} \left[ \mathrm{Li}_{2}(x) \left( \mathrm{Li}_{3}(x) - 4 \mathrm{Li}_{3}(x^{2}) \right) \right]_{0}^{1} + \frac{1}{8} \int_{0}^{1} \frac{\mathrm{Li}_{2}(x)}{x} \left( \mathrm{Li}_{2}(x) - 8 \mathrm{Li}_{2}(x^{2}) \right) \, dx \\
&= \frac{3}{8} \zeta(2)\zeta(3) + \frac{1}{8} \left[ \mathrm{Li}_{3}(x) \left( \mathrm{Li}_{2}(x) - 8 \mathrm{Li}_{2}(x^{2}) \right) \right]_{0}^{1} - \frac{1}{8} \int_{0}^{1} \frac{\mathrm{Li}_{3}(x)}{x} \left( \mathrm{Li}_{1}(x) - 16 \mathrm{Li}_{1}(x^{2}) \right) \, dx \\
&= -\frac{1}{2}\zeta(2)\zeta(3) + \frac{1}{8} \int_{0}^{1} \frac{\mathrm{Li}_{3}(x)}{x} \left( \log(1 - x) - 16 \log(1 - x^{2}) \right) \, dx \\
&= -\frac{1}{2}\zeta(2)\zeta(3) + \frac{1}{8} \int_{0}^{1} \frac{\mathrm{Li}_{3}(x) \log(1 - x)}{x} \, dx - 2 \int_{0}^{1} \frac{\mathrm{Li}_{3}(x) \log(1 - x^{2}) }{x}\, dx.
\end{align*} Now applying the substitution $x \mapsto x^{2}$ to the first integral, from the identity $\mathrm{Li}_{3}(x^{2}) = 4 \mathrm{Li}_{3}(x) + 4 \mathrm{Li}_{3}(-x)$ it follows that
\begin{align*}
\frac{1}{8} \int_{0}^{1} \frac{\mathrm{Li}_{3}(x) \log(1 - x)}{x} \, dx
&= \frac{1}{4} \int_{0}^{1} \frac{\mathrm{Li}_{3}(x^{2}) \log(1 - x^{2})}{x} \, dx \\
&= \int_{0}^{1} \frac{ \{ \mathrm{Li}_{3}(x) + \mathrm{Li}_{3}(-x) \} \log(1 - x^{2})}{x} \, dx.
\end{align*} This finally gives
\begin{align*}
I_{2} = -\frac{1}{2}\zeta(2)\zeta(3) - \int_{-1}^{1} \frac{\mathrm{Li}_{3}(x) \log(1 - x^{2}) }{x}\, dx.
\end{align*} Plugging this back to ( \text{(2)} ), we obtain
\begin{align*}
I_{1} = \frac{5}{4}\zeta(2)\zeta(3) - \frac{7}{2}\zeta(5) + I_{3}, \tag{3}
\end{align*} where ( I{3} ) is given by
\begin{align*}
I_{3} = - \int_{-1}^{1} \frac{\mathrm{Li}_{3}(x) \log(1 - x^{2})}{x} \, dx.
\end{align*}
Step 3. Finally, we evaluate $I_{3}$ into a closed form and thus complete the calculation. Since the integrand is holomorphic on the unit disc, we can shift up the contour of integration to the clockwise semicircular arc:
\begin{align*}
I_{3}
&= -i \int_{\pi}^{0} \mathrm{Li}_{3}(e^{i\theta}) \log(1 - e^{2i\theta}) \, d\theta
= - \frac{1}{i} \int_{0}^{\pi} \mathrm{Li}_{3}(e^{i\theta}) \log(1 - e^{2i\theta}) \, d\theta
\end{align*} Now expanding with Taylor series and taking advantage of the fact that $I_{3}$ is real, we have
\begin{align*}
I_{3}
&= \sum_{n=1}^{\infty} \sum_{m=1}^{\infty} \frac{1}{n^{3}} \frac{1}{m} \int_{0}^{\pi} \{ \cos (n\theta) \sin (2m\theta) + \sin(n\theta) \cos(2m\theta) \} \, d\theta \\
&= \sum_{n=1}^{\infty} \sum_{m=1}^{\infty} \frac{1}{n^{3}} \frac{1}{m} \int_{0}^{\pi} \cos (n\theta) \sin (2m\theta) \, d\theta \\
&\qquad + \sum_{n=1}^{\infty} \sum_{m=1}^{\infty} \frac{1}{n^{3}} \frac{1}{m} \int_{0}^{\pi} \sin(n\theta) \cos(2m\theta) \, d\theta \\
&=: S_{1} + S_{2}.
\end{align*} To evaluate $S_{1}$, we note the following simple identity
\begin{align*}
\int_{0}^{\pi} \cos (n\theta) \sin (2m\theta) \, d\theta
= \frac{1 + (-1)^{n-1}}{2} \int_{0}^{\pi} \cos (\tfrac{1}{2} n\theta) \sin (m\theta) \, d\theta.
\end{align*} Then we can write
\begin{align*}
S_{1} &= \sum_{n=0}^{\infty} \frac{1}{(2n+1)^{3}} \int_{0}^{\pi} \left\{ \sum_{m=1}^{\infty} \frac{\sin (m\theta)}{m} \right\} \cos((n+\tfrac{1}{2})\theta) \, d\theta, \\
S_{2} &= \sum_{m=1}^{\infty} \frac{1}{m} \int_{0}^{\pi} \left\{ \sum_{n=1}^{\infty} \frac{\sin(n\theta)}{n^{3}} \right\} \cos(2m\theta) \, d\theta.
\end{align*} But it is not hard to check that the following Fourier expansion holds:
\begin{align*}
\sum_{n=1}^{\infty} \frac{\sin n\theta}{n} = \Im \sum_{n=1}^{\infty} \frac{e^{i\theta}}{n} = - \Im \log(1 - e^{i\theta}) = \frac{\pi-\theta}{2}, \quad 0 < \theta < \pi,
\end{align*} \begin{align*}
\sum_{n=1}^{\infty} \frac{\sin n\theta}{n^{3}} = \frac{\theta ^3}{12}-\frac{\pi \theta ^2}{4}+\frac{\pi ^2 \theta }{6}, \quad 0 < \theta < \pi.
\end{align*} Plugging these to $S_{1}$, $S_{2}$ and performing some tedious calculation, we obtain
\begin{align*}
S_{1} = \sum_{n=0}^{\infty} \frac{2}{(2n+1)^{5}} = \frac{31}{16} \zeta(5)
\quad \text{and} \quad
S_{2} = -\frac{\pi^{2}}{16} \sum_{m=1}^{\infty} \frac{1}{m^{3}} = -\frac{3}{8}\zeta(2)\zeta(3).
\end{align*} Therefore we have
\begin{align*}
I = \frac{7}{8}\zeta(2)\zeta(3) - \frac{25}{16} \zeta(5).
\end{align*}
