Why the projection matrix can adjust the coordinates correctly?

If you target OpenGL ES 1.0 there is the OpenGL fixed-function pipeline. And here we have the GL_PROJECTION and the GL_MODELVIEW matrix stacks (actually there is another one for texture operations, but that is irrelevant here).

What happens is that the fixed-function pipeline simply concatenates/multiplies the top matrix of those two matrix stacks (let's refer to those as P and MV for projection matrix and modelview matrix, respectively) as P * MV during each draw call and transforms any drawn vertex v by P * MV * v.

When you target OpenGL ES 2.0 or higher, there are no such fixed-function matrix stacks anymore, but you would be using shaders and uploading your own matrices (e.g. as uniforms or by other means). There, you are completely free as to how you decomposition the final transformation (previously P * MV in the fixed-function pipeline). In the end, all that you want to do is to transform vertices from their space (typically called "model space") into "clip space" coordinates. Clip space is the only thing that OpenGL ES 2.0 (or higher) actually cares about.

It's just that decompositing the final matrix (model -> clip space) into projection, view and model matrices is convenient, because of the frequency that they change during a whole scene render. The projection matrix would typically only contain those transformations that stay constant and only needs to be refreshed when you resize the window (or reconfigure the vertical field of view). That is to say, this matrix transforms from "view space" to "clip space" and represents the intrinsic properties of a camera.

Now, in most cases there is also a dedicated "view" matrix, which represents the external properties of a camera, which is typically the position and orientation of the camera in the world. This matrix thusly changes when you rotate/move the virtual camera. And this matrix transforms from "model" to "view" space.

Next, there would be the model matrices, which represent the transformations of each model in the world.

And finally, when you want to render anything, you'd use the decomposition P * V * M (projection matrix multiplied by view matrix multiplied by model matrix) to transform a vertex through the whole chain of model -> clip space.

Now to the question about why we do the scaling in the projection matrix:

The reason is that once you are in view space (after the GL_MODELVIEW transformation or your custom V * M transformation), you have a space whose axes are aligned with the screen. So, you know what to scale (the X or Y axes) and by how much. You cannot simply scale the model, since that might also be rotated and translated arbitrarily.