You're correct that multiple linear layers can be equivalent to a single linear layer. As the other answers have said, a nonlinear activation function allows nonlinear classification. Saying that a classifier is nonlinear means that it has a nonlinear decision boundary. The decision boundary is a surface that separates the classes; the classifier will predict one class for all points on one side of the decision boundary, and another class for all points on the other side.
Let's consider a common situation: performing binary classification with a network containing multiple layers of nonlinear hidden units and an output unit with a sigmoidal activation function. $y$ gives the output, $h$ is a vector of activations for the last hidden layer, $w$ is a vector of their weights onto the output unit, and $b$ is the output unit's bias. The output is:
$$y = \sigma(hw + b)$$
where $\sigma$ is the logistic sigmoid function. Output is interpreted as the probability that the class is $1$. The predicted class $c$ is:
$$c = \left \{ \begin{array}{cl}
0 & y \le 0.5 \\
1 & y > 0.5 \\
\end{array} \right . $$
Let's consider the classification rule with respect to the hidden unit activations. We can see that the hidden unit activations are projected onto a line $hW + b$. The rule for assigning a class is a function of $y$, which is monotonically related to the projection along the line. The classification rule is therefore equivalent to determining whether the projection along the line is less than or greater than some threshold (in this case, the threshold is given by the negative of the bias). This means that the decision boundary is a hyperplane that's orthogonal to the line, and intersects the line at a point corresponding to that threshold.
I said earlier that the decision boundary is nonlinear, but a hyperplane is the very definition of a linear boundary. But, we've been considering the boundary as a function of the hidden units just before the output. The hidden unit activations are a nonlinear function of the original inputs, due to the previous hidden layers and their nonlinear activation functions. One way to think about the network is that it maps the data nonlinearly into some feature space. The coordinates in this space are given by the activations of the last hidden units. The network then performs linear classification in this space (logistic regression, in this case). We can also think about the decision boundary as a function of the original inputs. This function will be nonlinear, as a consequence of the nonlinear mapping from inputs to hidden unit activations.
This blog post shows some nice figures and animations of this process.