Using Taylor
\begin{align*}
f(X)=f(X_i)+f^{\prime}(X_i)(X-X_i)+\dfrac{f^{\prime \prime}(X_i)}{2}(X-X_i)^{2}+O((X-X_i)^{3})
\end{align*}
where $f(X)=X^3$, $X=B_{t_i}$ and $X_i=B_{t_{i-1}}$. Then,
\begin{align*}
B_{t_i}^{3}& =B_{t_{i-1}}^{3}+3B_{t_{i-1}}^{2}(B_{t_{i}}-B_{t_{i-1}})+3B_{t_{i-1}}(B_{t_{i}}-B_{t_{i-1}})^{2}+O((B_{t_{i}}-B_{t_{i-1}})^{3}) \\
%&=B_{t_{i-1}}^{3}+3B_{t_{i-1}}^{2}(B_{\Delta t_{i}})+3B_{t_{i-1}}(B_{\Delta t_{i}}^{2})+O(B_{\Delta t_{i}}^{3}) \\
&=B_{t_{i-1}}^{3}+3B_{t_{i-1}}^{2}(\Delta B_{t_{i}})+3B_{t_{i-1}}(\Delta B_{ t_{i}})^{2}+O((\Delta B_{t_{i}})^{3})
\end{align*}
Sum on $\Pi_n$ a partition of $[0,t]$, we have
\begin{align*}
\sum_{\Pi_n} B_{t_i}^{3}-B_{t_{i-1}}^{3}=\sum_{\Pi_n} 3B_{t_{i-1}}^{2}(\Delta B_{t_{i}})+3B_{t_{i-1}}(\Delta B_{ t_{i}})^{2}+O((\Delta B_{t_{i}})^{3}),
\end{align*}
adding and subtracting
\begin{align*}
\sum_{\Pi_n} B_{t_i}^{3}-B_{t_{i-1}}^{3}=\sum_{\Pi_n} 3B_{t_{i-1}}^{2}(\Delta B_{t_{i}})+3B_{t_{i-1}}(\Delta B_{ t_{i}})^{2} - 3B_{t_i}\Delta t_i + 3B_{t_i}\Delta t_i + O((\Delta B_{t_{i}})^{3}).
\end{align*}
So,
\begin{align*}
\lim_{n \to \infty} \sum_{\Pi_n} B_{t_i}^{3}-B_{t_{i-1}}^{3} = & \lim_{n \to \infty} \sum_{\Pi_n} 3B_{t_{i-1}}^{2}(\Delta B_{t_{i}})+ \lim_{n \to \infty} \sum_{\Pi_n}3B_{t_{i-1}}[(\Delta B_{ t_{i}})^{2} - \Delta t_i] \\
&+ \lim_{n \to \infty} \sum_{\Pi_n}3B_{t_i}\Delta t_i + \lim_{n \to \infty} \sum_{\Pi_n}O((\Delta B_{t_{i}})^{3}).
\end{align*}
Also, note that $$ \lim_{n \to \infty} \sum_{\Pi_n} 3 B_{t_{i-1}}[(\Delta B_{ t_{i}})^{2} - \Delta t_i] \displaystyle \longrightarrow 0 \ \mbox{ in } L^{2}$$
and $$ \lim_{n \to \infty} \sum_{\Pi_n}O((\Delta B_{t_{i}})^{3}) \longrightarrow 0 \ \mbox{in Probability.}$$
Therefore,
\begin{align*}
B_{t}^{3} = \int_0^t 3B_{s}^{2} \ d B_{s} +\int_0^t 3B_{s}\ ds \ \Longrightarrow \ \int_0^t B_{s}^{2} \ dB_s = \dfrac{B_t^3}{3}- \int_0^t B_{s}\ ds.
\end{align*}