Natural remanent magnetization diverts the direction of the total magnetization vector from the earth’s magnetic field direction. It is important to determine the magnetization direction because the processing and inversion of magnetic usually need the magnetization direction as a priori information. The conventional strategy for estimating the magnetization direction is to compute the crosscorrelation between the reduced-to-the-pole (RTP) field and one magnetic magnitude transform, which is insensitive to the magnetization direction. A new theory and method are proposed that the total magnetization direction is determined by computing the multiple correlation coefficients between the RTP field and the multiple direction-insensitive magnitudes transforms, such as the total magnitude anomaly and the normalized source strength. The linear regression between the RTP field and the direction-insensitive magnitude transforms is established. Then, the correlation is computed by the RTP field and the constructed regressor. The proposed method is tested on synthetic of single and complex models, respectively. Then, the method is applied to the field data of the Yeshan region (eastern China) and the Black Hill Norite (southern Australia). The proposed multiple correlation is an extension of the conventional crosscorrelation method that provides an accurate and robust way to estimate the magnetization direction from the total field anomaly.