This is an important question about sampling without replacement. I therefore offer a mathematical analysis, formulas, examples, approximations, and working code. Far more could be said, but to keep this post to a reasonable length I don't provide all the details or mention the implications.
Because $99$ is too large to permit simple illustrations, let's look at the same problem but with smaller numbers.
Analysis
To begin, consider sampling $k=2$ (instead of $k=5$) distinct values, uniformly, from the numbers $1, 2, 3, \ldots, n=5$ (instead of $n=99$). As we well know, the number of these samples--each of which has the same chance of being selected--is
$$\binom{n}{k} = \frac{n(n-1)\cdots(n-k+2)(n-k+1)}{k(k-1)\cdots(2)(1)},$$
which is equal to $10$ in this example. Here is an illustration of these ten possible samples.
Read these by counting the filled squares in each row, from bottom to top. For instance, the top left panel has one gray square in its bottom row and two gray squares in the next row above it. This indicates the sample $\{2,1\}.$
In the title of each panel I have indicated the total number of colored squares in the rectangles followed by their counts in the rows of the rectangles (from top to bottom).
Notice the following:
The gray triangle of squares at the left of each panel guarantees there will be at least one square in the bottom row and two squares in the top row.
The colored squares fill in the $2\times 3 = k \times n-k$ rectangle shown in each panel.
The number of colored squares in any row of the rectangle always equals or exceeds the number in any rows below it.
These ways of filling a rectangle are called Young tableaux. I call the ones in the figure "augmented" tableaux to indicate the padding by gray squares at the left (which is not standard): $k=2$ squares at the top, $k-1$ squares in the next row, and so on down to the bottom row of $1$ square. The total number of padding squares therefore is $k(k+1)/2.$
- Each panel thereby depicts a distinct subset of $n=5$ numbers. Why? Because
- The counts of colored squares in the rectangles do not decrease as you move upwards, yet
- The counts of gray squares always increase by $1,$ whence the total number of squares in the rows is strictly increasing.
- The largest number of squares in any row amounts to $k$ (gray squares in the top row) plus $n-k$ (the width of the rectangles), totaling $n.$
- There is at least one square in each row.
That's why there are $\binom{n}{k} = \binom{5}{2}=10$ panels in the figure.
- Different colors indicate different sums of the values in the subsets.
- The sums range from $k(k+1)/2$ at the upper left, counting all the gray squares, to $k(k+1)/2 + (n-k)k$ at the lower right. These correspond to the extreme sums $1+2+\cdots + k$ on the one hand and $(n-k+1) + (n-k+2) + \cdots + n$ on the other.
To help you get these ideas intuitively, here's a plot of all $\binom{6}{3}=20$ augmented tableaux for samples of size $k=3$ from the numbers $1, 2, \ldots, n=6.$
Using this figure you can read off the probability distribution of the sum of three fair, ordinary dice conditional on all three values being different: just add $1+2+3=6$ to each sum shown.
Because each augmented tableau has an equal chance of $\binom{n}{k}$ of being selected,
The chance that the sum of values in a sample of size $k$ taken from $1,2,\ldots,n$ without replacement is equal to $r$ is $1/\binom{n}{k}$ times the count of Young tableaux of dimension $k\times n-k$ containing $r - k(k+1)/2$ squares.
It is convenient to roll all these counts up into a polynomial. Its powers are the sums and its coefficients are the numbers of samples in which each sum occurs. For example, the figure shows one sum of $0,$ one of $1,$ two sums (green) of $2,$ and so on. The polynomial is
$$\binom{5}{2}_q = 1 q^0 + 1 q^1 + 2 q^2 + 2 q^3 + 2 q^4 + 1 q^5 + 1 q^6.$$
It is known as the Gaussian binomial coefficient. As you can see in the figure,
The Gaussian binomial coefficient $\binom{n}{k}_q$ shows how the samples of size $k$ from the numbers $1,2,\ldots, n$ break down by their sums relative to the smallest possible sum of $k(k+1)/2.$
Another way to put this is that
you can read off the probability the sum will equal $r$ by inspecting the coefficient of $q^r$ in the polynomial $$q^{k(k+1)/2}\frac{\binom{n}{k}_q}{\binom{n}{k}}.$$
This means $q^{k(k+1)/2}\binom{n}{k}_q/\binom{n}{k}$ is the probability generating function (pgf) for the distribution of the sums.
Computation
A great deal is known about Young tableaux: they are important in combinatorics. One of the basic ideas is that you can generate them recursively. If the bottom row of a $k\times n$ tableau is empty, it's just a $k-1\times n$ tableau with an extra empty row. Otherwise, all its rows are filled, which means the left column of the rectangle is entirely filled. When you lop off that column you are left with a $k\times n-1$ tableau. Thus, by adjoining an empty bottom row or full left column as appropriate, you can find all the $k\times n$ tableaux once you have found all the $k-1\times n$ and $k\times n-1$ tableaux.
This yields a recurrence relation for the probability generating functions as well as formulas for them. In particular,
$$\binom{n}{k}_q = \frac{[n]_q [n-1]_q \cdots [n-k+1]_q}{[k]_q [k-1]_q \cdots [1]_q}$$
where the "$q$-numbers" are given by
$$[k]_q = \frac{1 - q^k}{1-q} = 1 + q + q^2 + \cdots + q^{k-1}.$$
(Compare this to the formula I gave earlier for $\binom{n}{k}:$ it is identical, but $q$-numbers are used rather than ordinary integers.)
Continuing our example,
$$\begin{aligned} \binom{5}{2}_q &= \frac{[5]_q [4]_q}{[2]_q [1]_q} = \frac{(1 + q + q^2 + q^3 + q^4)(1+q+q^2+q^3)}{(1+q)(1)} \\ &= 1 + q + 2q^2 + 2q^3 + 2q^4 + q^5 + q^6.\end{aligned}$$
Adding $k(k+1)/2 = 3$ to each exponent and dividing by $\binom{5}{2}=10$ gives the pgf
$$\frac{1}{10}\left(q^3 + q^4 + 2q^5 + 2q^6 + 2q^7 + q^8 + q^9\right).$$
For instance, the chance the sum will equal $7$ is the coefficient of $q^7,$ equal to $2/10 = 1/5.$
This leads to several possible algorithms for computing the pgf. You can apply the recursion directly (just like computing Pascal's Triangle, but the entries are now polynomials); or you can apply the formula, using either convolution and deconvolution via the Fast Fourier Transform or direct implementation of polynomial multiplication and division.
I tried all these variations. Convolution is fastest (as one would expect), but it runs into severe numerical problems by the time $n=99$ due to instability in the division (deconvolution). By limiting the FFT to multiplication and performing the division directly you can avoid most of these problems, but then the algorithm isn't appreciably faster than the recursion. I have therefore used the recursion in the following calculations.
Comparing polynomial multiplication to addition of random variables leads to a simple Normal approximation. The expectation of the sum is
$$\mu(n,k) = k\frac{n+1}{2},$$
of course. Its variance is
$$\sigma^2(n,k) = \frac{k}{n} \left[\frac{n-k}{n-1} \sum_{i=1}^n i^2 - \frac{n-k}{n(n-1)} \left(\sum_{i=1}^n i\right)^2\right].$$
(Although the sums are readily evaluated in closed form, I have left them there to indicate how this formula was derived.)
The probability assigned by the Normal distribution with parameters $(\mu,\sigma)$ to the interval $[r-1/2,r+1/2]$ approximates the chance the sum will equal $r.$
Application
Let's answer the question.
The recursive method computes all the Gaussian binomial coefficients for $n.$ Thus, having computed them for $n=99$ and $k=5,$ I also looked at the sum distributions for other sample sizes $k = 2, 3, 10.$ We would expect the distribution to become closer to Normal as $k$ increases.
Here are plots of the pgfs for these four values of $k$ (including the $k=5$ of the question) on which I have overplotted the Normal curve in red. ("WOR" in the titles means "without replacement"). The entire calculation took 0.01 seconds.
Clearly the Normal approximation is excellent by the time $k$ is as large as $10.$ But because we often are interested in extreme probabilities and they can scarcely be seen in the top row of plots, I have also plotted the differences in quantiles between the true sum distributions and their Normal approximations in the lower strip. They all use the same scale so you can see how the errors decrease as $k$ grows.
Code
I'm sure someone (you know who you are!) will request the code I used to do these calculations. qBinomial
, reproduced below, implements the recursive method. Its output is the coefficients of the Gaussian binomial coefficient $\binom{n}{k}_q$ in order of increasing powers from $0$ on up. If you do not specify $k,$ it returns a list of all such coefficients for $k=0, 1, \ldots, n.$
Since there is an (obvious) symmetry in the coefficients arising from the fact that taking a sample of size $k$ of sum $r$ leaves a complementary sample of size $n-k$ of sum $n(n+1)/2-r,$ the list is indexed by the pairs $\{k,n-k\}.$ This is revealed in the names in this output for $n=5:$
> qBinomial(5)
$`0 or 5`
[1] 1
$`1 or 4`
[1] 1 1 1 1 1
$`2 or 3`
[1] 1 1 2 2 2 1 1
That last entry summarizes the first figure.
If you would like individual probabilities, normalize the polynomial and index into it with an offset of $k(k+1)/2.$ Here, for instance, is the chance that the sum equals $220:$
> p <- qBinomial(99, 5) / choose(99, 5)
> p[[220 - 5*(5+1)/2]]
[1] 0.00554545
(You can inspect the output of qBinomial
before normalization to find the exact probability is $396628/71523144.$)
Here is the chance the sum equals or exceeds $400:$
> sum(p[seq(400 - 5*(5+1)/2, length(p))])
[1] 0.007088153
This is the code. Running time is theoretically cubic in $n$ but in R
it seems to scale like $O(n^5).$ It works well enough to compute values for $n \le 200$ quickly.
#
# The Gaussian q-Binomial coefficients, returned as an array indexed (from 0)
# for the powers of q. When `k` is absent, it will return all the coefficients
# for `n` as a list # of arrays covering the first half. Use the k <--> n-k
# symmetry to index into it (see the last line for an example).
#
# Example:
# n <- 6
# X.lst <- qBinomial(n)
# k <- 4
# (X <- X.lst[[min(k, n-k) + 1]])
#
# Output:
# [1] 1 1 2 2 3 2 2 1 1
#
qBinomial <- function(n, k) {
if (n==0) return(list(`0`=c(1)))
#
# Shift `x` to the right by `r` and add to `y`.
#
addWithShift <- function(x, y, r) {
n <- max(length(x) + r, length(y))
c(y, rep(0, n-length(y))) + c(rep(0, r), x)
}
#
# Recursion.
#
R <- list(c(1), c(1)) # Values for r=0 through ceiling((m+1)/2)
if (n > 1) {
for (m in seq(2, n)) {
m2 <- ceiling((m+1)/2)
R. <- R
for (r in seq(1, m2-1)) R[[r+1]] <- addWithShift(R.[[r+1]], R.[[r]], r)
if(2*m2 == m+1) R[[m2+1]] <- R[[m2]]
}
}
#
# Document the indexing.
#
names(R) <- paste(seq(0, length(R)-1), n - seq(0, length(R)-1), sep=" or ")
if (missing(k)) return(R) else return(R[[1 + min(k, n-k)]])
}