You might find this short note helpful ($\LaTeX$ code available [1]
if the binary link breaks). I am reproducing the relevant part below:
Theorem. (Folklore) Learning an unknown distribution over a known domain of size $n$, up to total variation $\varepsilon\in(0,1]$, and with error probability $\delta\in(0,1]$, has sample complexity $O\!\left(\frac{n+\log(1/\delta)}{\varepsilon^2}\right)$. (Moreover, this can be done efficiently.)
Proof. Consider the empirical distribution $\tilde{p}$ obtained by drawing $m$ independent samples $s_1,\dots,s_m$ from the underlying distribution $p\in\Delta([n])$:
\begin{equation}\label{def:empirical}
\tilde{p}(i) = \frac{1}{m} \sum_{j=1}^m \mathbb{1}_{\{s_j=i\}}, \qquad i\in [n]
\end{equation}
First, we bound the expected total variation distance between $\tilde{p}$ and $p$, by using $\ell_2$ distance as a proxy:
$$
\mathbb{E}{ d_{\rm TV}(p,\tilde{p}) }
=\frac{1}{2}\mathbb{E}{ \lVert{p-\tilde{p}}\rVert_1}
=\frac{1}{2}\sum_{i=1}^n\mathbb{E}{ \lvert{p(i)-\tilde{p}(i)}\rvert }
\leq\frac{1}{2}\sum_{i=1}^n\sqrt{\mathbb{E}{ (p(i)-\tilde{p}(i))^2} }
$$
the last inequality by Jensen. But since, for every $i\in[n]$, $m\tilde{p}(i)$ follows a $\operatorname{Bin}({m},{p(i)})$ distribution, we have
$\mathbb{E}{ (p(i)-\tilde{p}(i))^2} = \frac{1}{m^2}\operatorname{Var}[m\tilde{p}(i)] = \frac{1}{m}p(i)(1-p(i))$, from which
$$
\mathbb{E}{ d_{\rm TV}(p,\tilde{p}) } \leq\frac{1}{2\sqrt{m}}\sum_{i=1}^n\sqrt{p(i)} \leq \frac{1}{2}\sqrt{\frac{n}{m}}
$$
the last inequality this time by Cauchy—Schwarz. Therefore, for $m\geq \frac{n}{\varepsilon^2}$ we have $\mathbb{E}{ d_{\rm TV}(p,\tilde{p}) }\leq \frac{\varepsilon}{2}$.
Next, to convert this expected result to a high probability guarantee, we apply McDiarmid's inequality to the random variable $f(s_1,\dots,s_m) \stackrel{\rm def}{=} d_{\rm TV}(p,\tilde{p})$, noting that changing any single sample cannot change its value by more than $c\stackrel{\rm def}{=} 1/m$:
$$
\mathbb{P}\left\{ \lvert{f(s_1,\dots,s_m) - \mathbb{E}{f(s_1,\dots,s_m)}\rvert} \geq \frac{\varepsilon}{2} \right\} \leq 2e^{-\frac{2\left(\frac{\varepsilon}{2}\right)^2}{mc^2}} = 2e^{-\frac{1}{2}m\varepsilon^2}
$$
and therefore as long as $m\geq \frac{2}{\varepsilon^2}\ln\frac{2}{\delta}$, we have $\lvert{f(s_1,\dots,s_m) - \mathbb{E}{f(s_1,\dots,s_m)}\rvert} \leq \frac{\varepsilon}{2}$ with probability at least $1-\delta$. $\square$
There is a second proof, somewhat more fun, given in that short note (credit to John Wright for pointing it out, and emphasizing it's the "fun" one). Here it is:
Proof. Again, we will analyze the behavior of the empirical distribution $\tilde{p}$ over $m$ i.i.d. samples from the unknown $p$. Recalling the definition of total variation distance, note that $d_{\rm TV}({p,\tilde{p}}) > \varepsilon$ literally means there exists a subset $S\subseteq [n]$ such that $\tilde{p}(S) > p(S) + \varepsilon$. There are $2^n$ such subsets, so we can do a union bound.
Fix any $S\subseteq[n]$. We have
$$
\tilde{p}(S) = \tilde{p}(i) = \frac{1}{m} \sum_{i\in S} \sum_{j=1}^m \mathbb{1}_{\{s_j=i\}}
$$
and so, letting $X_j \stackrel{\rm def}{=} \sum_{i\in S}\mathbb{1}_{\{s_j=i\}}$ for $j\in [m]$, we have
$
\tilde{p}(S) = \frac{1}{m}\sum_{j=1}^m X_j
$ where the $X_j$'s are i.i.d. Bernoulli random variable with parameter $p(S)$. Then, by a Chernoff bound (actually, Hoeffding):
$$
\mathbb{P}\left\{ \tilde{p}(S) > p(S) + \varepsilon \right\} = \mathbb{P}\left\{ \frac{1}{m}\sum_{j=1}^m X_j > \mathbb{E}\left[\frac{1}{m}\sum_{j=1}^m X_j\right] + \varepsilon \right\} \leq e^{-2\varepsilon^2 m}
$$
and therefore $\mathbb{P}\left\{ \tilde{p}(S) > p(S) + \varepsilon \right\} \leq \frac{\delta}{2^n}$ for any $m\geq \frac{n\ln 2+\log(1/\delta)}{2\varepsilon^2}$. A union bound over these $2^n$ possible sets $S$ concludes the proof:
$$
\mathbb{P}\left\{ \exists S\subseteq [n] \text{ s.t. }\tilde{p}(S) > p(S) + \varepsilon \right\} \leq 2^n\cdot \frac{\delta}{2^n} = \delta
$$
and we are done. $\square$
Note: a lower bound of $\Omega(\frac{n}{\varepsilon^2})$ (also folklore) is easy to derive from Assouad's lemma, by considering the family of distributions over $[n]$ where each pair of consecutive elements $(2i,2i+1)$ has either probabilities $(\frac{1+c\varepsilon}{n},\frac{1-c\varepsilon}{n})$ or $(\frac{1-c\varepsilon}{n},\frac{1+c\varepsilon}{n})$ for some suitable constant $c>0$. (Intuitively and a bit misleadingly: any learning algorithm has to "figure out" at least $\Omega(n)$ of these independent choices, but each of them requires $\Omega(1/\varepsilon^2)$ samples.)
[1] Public GitHub: https://github.com/ccanonne/probabilitydistributiontoolbox (includes the source for the note on Assouad's lemma as well).
