Intuitively understand why the Poisson distribution is the limiting case of the binomial distribution In "Data Analysis" by D. S. Sivia, there is a derivation of the Poisson distribution, from the binomial distribution. 
They argue that the Poisson distribution is the limiting case of the binomial distribution when $M\rightarrow\infty$, where $M$ is the number of trials.
Question 1: How can that argument intuitively be understood?
Question 2: Why is the large-$M$ limit of $\frac{M!}{N!(M-N)!}$ equal to $\frac{M^{N}}{N!}$, where $N$ is the number of successes in $M$ trials? (This step is used in the derivation.)
 A: Let me provide an alternate heuristic.  I'm going to show how to approximate the Poisson process as a binomial (and argue that the approximation is better for many trials with low probability).  Therefore the binomial distribution must tend to the Poisson distribution.
Let's say events are happening with a constant rate in time.  We want to know the distribution of how many events happened in a day, knowing that the expected number of events is $\lambda$.
Well, the expected number of events per hour is $\lambda/24$.  Let's pretend that this means that the probability of an event happening in a given hour is $\lambda/24$.  [it's not quite right, but it is a decent approximation if $\lambda/24 \ll 1$ basically if we can assume multiple events don't happen in the same hour].  Then we can approximate the distribution of the number of events as a binomial with $M=24$ trials, each having success probability $\lambda/24$.
We improve the approximation by switching our interval to minutes.  Then it's $p=\lambda/1440$ with $M=1440$ trials.  If $\lambda$ is around, say 10, then we can be pretty confident that no minute had two events.
Of course it gets better if we switch to seconds.  Now we're looking at $M=86400$ events each with the small probability $\lambda/86400$.  
No matter how big your $\lambda$ is, I can eventually choose a small enough $\Delta t$ such that it's very likely that no two events happen in the same interval.  Then the binomial distribution corresponding to that $\Delta t$ will be an excellent match to the true Poisson distribution.  
The only reason they aren't exactly the same is that there is a non-zero probability that two events happen in the same time interval.  But given there are only around $\lambda$ events and they are distributed into some number of bins much greater than $\lambda$, it's unlikely that any two of them lie in the same bin.
Or in other words, the binomial distribution tends to the Poisson distribution as $M \to \infty$ if the success probability is $p=\lambda/M$.
A: The problem is that your characterization of the Poisson as a limiting case of the binomial distribution is not quite correct as stated.
The Poisson is a limiting case of the binomial when:  $$M \to \infty \quad \color{red}{\text{and} \quad Mp \to \lambda.}$$  The second part is important.  If $p$ remains fixed, the first condition implies that the rate will also increase without bound.
What the Poisson distribution assumes is that events are rare.  What we mean by "rare" is not that the rate of events is small--indeed, a Poisson process may have a very high intensity $\lambda$--but rather, that the probability of an event occurring at any instant in time $[t, t + dt)$ is vanishingly small.  This is in contrast to a binomial model where the probability $p$ of an event (e.g. "success") is fixed for any given trial.
To illustrate, suppose we model a series of $M$ independent Bernoulli trials each with probability of success $p$, and we look at what happens to the distribution of the number of successes $X$ as $M \to \infty$.  For any $N$ as large as we please, and no matter how small $p$ is, the expected number of successes $\operatorname{E}[X] = Mp > N$ for $M > N/p$.  Put another way, no matter how unlikely the probability of success, eventually you can achieve an average number of successes as large as you please if you perform sufficiently many trials.  So, $M \to \infty$ (or, just saying "$M$ is large") is not enough to justify a Poisson model for $X$.
It is not difficult to algebraically establish $$\Pr[X = x] = e^{-\lambda} \frac{\lambda^x}{x!}, \quad x = 0, 1, 2, \ldots$$ as a limiting case of $$\Pr[X = x] = \binom{M}{x} p^x (1-p)^{M-x}, \quad x = 0, 1, 2, \ldots, M$$ by setting $p = \lambda/M$ and letting $M \to \infty$.  Other answers here have addressed the intuition behind this relationship and provided computational guidance as well.  But it is important that $p = \lambda/M$.  You can't ignore this.
A: Question 1
Recall the definition of the binomial distribution:

a frequency distribution of the possible number of successful outcomes in a given number of trials in each of which there is the same probability of success.

Compare this to the definition of the Poisson distribution:

a discrete frequency distribution which gives the probability of a number of independent events occurring in a fixed time.

The substantial difference between the 2 is the binomial is in $n$ trials, Poisson is over a time period $t$. How can the limit occur intuitively?
Lets say that you have to keep running Bernoulli trials for all eternity. Moreover, you run $n = 30$ per minute. Per minute you count each success. So for all eternity you are running a $Bin(p,30)$ process every minute. Over 24 hours, you have a $Bin(p,43200)$.
As you get tired, you are asked "how many successes occurred between 18:00 and 19:00?". You're answer might be $30*60*p$, i.e. you provide the average successes in an hour. That sounds a lot like the Poisson parameter $\lambda$ to me.
A: Question 2) 
$$\frac{\frac{M!}{N!(M-N)!}}{\frac{M^N}{N!}} = \frac{M(M-1)\dots(M - N + 1)}{M^N} = 1(1 - \frac{1}{M})\dots(1 - \frac{N - 1}{M})$$
So taking the limit for fixed $N$
$$\lim_{M \to \infty} \frac{\frac{M!}{N!(M-N)!}}{\frac{M^N}{N!}} = \lim_{M \to \infty} 1(1 - \frac{1}{M})\dots(1 - \frac{N - 1}{M}) = 1$$
A: I will try a simple intuitive explanation.  Record that for a binomial random variable $X \sim \text{Bin}(n,p)$ we have expectation is $n p$ and variance is $n p (1-p)$.   Now think that $X$ records the number of events in a very large number $n$ of trials, each with a very small probability $p$, such that we are very close to $1-p=1$ (really $\approx$).  Then we have $np=\lambda$ say, and $n p (1-p) \approx n p 1 =\lambda$, so the mean and variance are both equal to $\lambda$.  Then remember that for a poisson distributed random variable, we always have mean and variance equal!  That is at least a plausibility argument for the poisson approximation, but not a proof.
Then look at it from another viewpoint, the poisson point process https://en.wikipedia.org/wiki/Poisson_point_process   on the real line.  This is the distribution of random points on the line that we gets if random points occur according to the rules:


*

*points in disjoint intervals are independent

*probability of a random point in a very short interval is proportional to length of interval

*probability of two or more points in a very short interval is essentially zero. 


Then the distribution of number of points in a given interval (not necessarily short) is Poisson (with parameter $\lambda$ proportional to length).  Now, if we divide this interval in very many, equally very short subintervals ($n$), the probability of two or more points in a given subinterval is essentially zero, so that number will have, to a very good approximation, a bernolli distribution, that is, $\text{Bin}(1,p)$, so the sum of all this will be $\text{Bin}(n,p)$, so a good approximation of the poisson distribution of number of points in that (long) interval.
Edit from @Ytsen de Boer (OP): question number 2 is satisfactorily answered by @Łukasz Grad.
A: I can only attempt a part answer and it is about the intuition for Question 2, not a rigorous proof.
The binomial coefficient gives you the number of samples of size $N$, from $M$, without replacement and without order.
Here though $M$ becomes so large that you may approximate the scenario as sampling with replacement in which case you get
$M^N$ ordered samples. If you don't care about the order of the $N$ objects chosen this reduces to $M^N/N!$ because those $N$ objects 
can be ordered in $N!$ ways.
A: 
I think this is the best example that intuitively explains how binomial distribution converges to normal with large number of balls. Here, each ball has equal probability of falling on either side of the peg in each layer and all the balls have to face same number of pegs. It can be easily seen that as the number of balls goes very high the distribution of balls in different sections will be like normal distribution.
My answer to your question 2 is same as the answer given by Lukasz.
