The method of moments estimator is $\hat \theta_n = 2\bar X_n,$ and it is unbiased.
It has a finite variance (which decreases with increasing $n$)
and so it is also consistent; that is, it converges in probability to $\theta.$
I have not checked your proof of consistency, which seems inelegant and incorrect (for one thing, the $\epsilon$ disappears in the second line).
You should be able to use a straightforward application of Chebyshev's inequality to show that
$\lim_{n \rightarrow \infty}P(|\hat \theta_n - \theta| <\epsilon) = 1.$
However, $\hat \theta_n$ does not have the minimum variance among unbiased
estimators. The maximum likelihood estimator is the maximum of the $n$
values $X_i$ (often denoted $X_{(n)}).$ The estimator $T = cX_{(n)},$ where $c$ is constant depending on $n,$ is unbiased and has minimum variance among
unbiased estimators (UMVUE).
Both estimators are illustrated below for $n = 10$ and $\theta = 5$ by
simulations in R statistical software. With a 100,000 iterations
means and variances should be accurate to about two places. They
are not difficult to find analytically.
m = 10^5; n = 10; th = 5
x = runif(m*n, 0, th)
DTA = matrix(x, nrow=m) # m x n matrix, each row a sample of 10
a = rowMeans(DTA) # vector of m sample means
w = apply(DTA, 1, max) # vector of m maximums
MM = 2*a; UMVUE = ((n+1)/n)*w
mean(MM); var(MM)
## 5.003658 # consistent with unbiasedness of MM
## 0.8341769 # relatively large variance
mean(UMVUE); var(UMVUE)
## 5.002337 # consistent with unbiasedness of UMVUE
## 0.207824 # relatively small variance
The histograms below illustrate the larger variance of the method of
moments estimator.
Your calculation of $E(Y_n)$ is wrong, at the very last calculus step. You can check that $P(Y_n<X) = 1$, so $E(Y_n)=X$ is immediately suspect. (By the "Lake Woebegone" principle.) The correct integral is $E(Y_n)=n/(n+1)$.
That mistake aside, your general plan for figuring out the density function for $Y_n$ and expressing its moment with an integral is sound. You should be able to get the variance of $Y_n$ the same way.
Best Answer
Consider $$X_{(n)}=\max(X_1,X_2,...,X_n)$$ We have $$P(X_{(n)}< m)=P(\cap_{k\leq n}\{X_k<m\})=(1-p^m)^n\implies \sum_{n \in \mathbb{N}}{\underbrace{(1-p^m)}_{<1}}^n<\infty$$ By Borel-Cantelli I we have $P(X_{(n)}<m\textrm{ i.o.})=0$. This means that for a.a. $\omega,\exists N(\omega):\,\forall n \geq N(\omega)$ we have $X_{(n)}(\omega)=m$. But this means $X_{(n)}$ converges a.a. to $m$: a.a. convergence implies convergence in probability, so the estimator is consistent.