Why Vgs-Vth is not equal to vdsat?

First of all, MOSFET model used today is pretty complicated. They can not be fully-expressed in simple equation such as
\( I_{D} = \dfrac{1}{2}\mu C_{ox} \dfrac{W}{L} \left( V_{GS} - V_{TH}\right)^2\left( 1+ \lambda V_{DS}\right) \)
Therefore, \( V_{DSAT} \) is not exactly same as \( V_{GS} - V_{TH} \). Try to google something like BSIM4.
+ Due to sub-threshold behvaior of MOSFET, notion of \( V_{DSAT} \) itself is also very ambiguous even in the simple equation.

I want to use the current formula in the saturation region, which of the following is correct?

Click to expand...

Depending on the channel size, there will be difference on the formula; in short; you need to go a little further in theory to get the exact answer.

First of all, drain saturation voltage and overdrive voltage are two completely different quantities.
Based on some physical model (can be ACM, can be EKV) derived from relationship of channel's charge profile along the transistor we can clearly get "ideal" curves bounding Vod and VDSsat versus drain current (or inversion coefficient as somebody prefer). Below is such a plot for ideal transistor (no short channel effects, ideally flat output characteristic) at room temperature:

There is only one single point for which both quantities are equal.
Here, VDSsat has been defined as point for which forward current is 100× higher then reverse one. This is of course the matter of convention, however, despite of this choice both relationship (for Vod and VDSSat) are always different.