Abnormal value of the multibeam echosounder system is the main influencing factor of accuracy of hydrographic survey and mapping accuracy of seabed topography. Good functional and stochastic models for multibeam bathymetric data processing are necessary to protect accurate display of multibeam bathymetric data in seabed and safety navigation of ships. This study proposed a robust least square collocation method based on the improved multi-quadric function to represent the variation trend of seabed topography accurately and disclose the real seabed topography along with the systematic and random problems of nonsteady multibeam bathymetric data of seabed topography. Fitting methods of seabed topography trend surface and the elimination of trend items were analyzed. The solving method of the empirical covariance function in the least square collocation and the method to eliminate abnormal values in multibeam bathymetric data were proposed. Validity of the proposed method was proved by an experiment. Results demonstrate that the fitted trend terms in the improved multi-quadric function can reflect the overall trend changes of the study area effectively. The applied robust least square collocation method can effectively reflect the systematicness and randomness of nonsteady bathymetric data of seabed topography. This method can overcome the influences of gross errors or abnormal points on the multibeam bathymetric data and can thus increase the accuracy of the estimation points. Compared with the traditional quadratic polynomial surface fitting method, the proposed method is superior in reflecting the overall trend changes in seabed topography, quality of multibeam bathymetric data, and mapping accuracy of seabed topography. This study can increase the detection capacity and effects of seabed topography at a high working efficiency and lay foundations for scientific research and marine exploration in the future.