Assume that $E[( \theta-\hat{\theta}_n)^2] = O(n^{-1/2})$. How can I show that $\frac{1}{n}\sum_{i=1}^n (\theta-\hat{\theta}_n)^2 = O_p(n^{-1/2})?$
What I'm trying to ask is: if the expected value of some quantity is $O(n^{-1/2})$, does that imply the sample mean of that quantity is also $O_p(n^{-1/2})$?
You want to show that for any $\varepsilon>0$, there exist two positive reals $M$ and $N$ such that $$\mathbb P\left(\frac{1}{\sqrt n}\sum_{i=1}^n (\theta-\hat{\theta}_n)^2 > M\right)\le\varepsilon\ \ \forall n\ge N $$ Let $\varepsilon>0$. For any $n\ge1$ and $M>0$ you have, by Markov inequality, that $$\begin{align}\mathbb P\left(\frac{1}{\sqrt n}\sum_{i=1}^n (\theta-\hat{\theta}_n)^2 > M\right)&\le\frac 1 M\mathbb E\left[\frac{1}{\sqrt n}\sum_{i=1}^n (\theta-\hat{\theta}_n)^2 \right]\\ &\le\frac{1}{M\sqrt n}\mathbb E\left[\sum_{i=1}^n (\theta-\hat{\theta}_n)^2 \right]\\ &\le\frac{\sqrt n}{M}\mathbb E\left[(\theta-\hat{\theta}_n)^2 \right]\\\end{align} $$ But, because you assume that $E[( \theta-\hat{\theta}_n)^2] = O(n^{-1/2})$ (it's a deterministic sequence so I assume the big O is the usual one, otherwise it makes no sense), you know that for any positive $M'$, there exists $N'$ such that, for all $n\ge N'$ : $$\sqrt nE[( \theta-\hat{\theta}_n)^2] \le M' $$ Now take $M'=M\varepsilon$ and you can conclude that the desired inequality holds for all $n\ge N'$.
[Note that I assumed that the $\theta$ in the sum were supposed to be $n$ i.i.d. realizations of $\theta$ called $\theta_i$.]
We can interpret the $O_p$ in terms of the quantile functions.
If
$$X_n = {O}_p(a_n) \quad \text{for }n \to \infty$$
then in terms of the quantile functions for $X_n$
$$\forall p \in (0,1): \lbrace Q_{X_n}(p) = {O}(a_n) \quad \text{for }n \to \infty \rbrace$$
See this also in the question: What is difference between $\hat{X}_n \overset{p}{\to} \bar{x}$ and $(\hat{X}_n - \bar{x}) = o_p(1)$?
If the quantile function is continuous, then we can express the expectation value as an integral over the quantile function
$$E[X_n] = \int_0^1 Q_{X_n}(p) \,\text{d}p$$
And for non-negative variables $X_n$ we have for every $p_c<1$
$$E[X_n] \geq \int_{p_c}^1 Q_{X_n}(p) \,\text{d}p \geq \int_{p_c}^1 Q_{X_n}(p_c) \,\text{d}p = (1-p_c) Q_{X_n}(p_c)$$
From which it follows that $Q_{X_n}(p_c) = O\left(E[X_n]\right)$ and the quantiles $p_c<1$ have at least the same order of convergence as the expectation value.
Say the quantiles function is $$Q_{X_n}(p) = \log(n) n^{p/2-1}$$ then $$E[X_n] = 2n^{-0.5} + 2n^{-n} = O\left(n^{-0.5}\right)$$
The quantile functions have convergence $$O\left(\log(n)n^{-a}\right)$$ with $$a = 0.5 \cdot (2-p)$$ For $p =1$ this power is $0.5$ and the convergence bound is $O\left(\log(n)n^{-0.5}\right)$, which means less fast convergence. However, for any $p<1$ the power is $a>0.5$ and the convergence is faster than the convergence of the expectation. We can only make examples where the quantiles $p=1$ converges less fast.
