diff --git a/doc/api/index.rst b/doc/api/index.rst
index e08f1c0aacb..4d9a2616d96 100644
--- a/doc/api/index.rst
+++ b/doc/api/index.rst
@@ -167,12 +167,13 @@ and store them in the GMT cache folder.
     :toctree: generated
 
     datasets.load_earth_relief
+    datasets.load_fractures_compilation
+    datasets.load_hotspots
     datasets.load_japan_quakes
+    datasets.load_mars_shape
     datasets.load_ocean_ridge_points
     datasets.load_sample_bathymetry
     datasets.load_usgs_quakes
-    datasets.load_fractures_compilation
-    datasets.load_hotspots
 
 .. automodule:: pygmt.exceptions
 
diff --git a/pygmt/datasets/__init__.py b/pygmt/datasets/__init__.py
index 1ac58d23b87..141e7238804 100644
--- a/pygmt/datasets/__init__.py
+++ b/pygmt/datasets/__init__.py
@@ -7,6 +7,7 @@
     load_fractures_compilation,
     load_hotspots,
     load_japan_quakes,
+    load_mars_shape,
     load_ocean_ridge_points,
     load_sample_bathymetry,
     load_usgs_quakes,
diff --git a/pygmt/datasets/samples.py b/pygmt/datasets/samples.py
index 18ce24f0a30..daf89014f06 100644
--- a/pygmt/datasets/samples.py
+++ b/pygmt/datasets/samples.py
@@ -149,3 +149,28 @@ def load_hotspots():
     columns = ["longitude", "latitude", "symbol_size", "place_name"]
     data = pd.read_table(filepath_or_buffer=fname, sep="\t", skiprows=3, names=columns)
     return data
+
+
+def load_mars_shape():
+    """
+    Load a table of data for the shape of Mars.
+
+    This is the ``@mars370d.txt`` dataset used in GMT examples, with data and
+    information from Smith, D. E., and M. T. Zuber (1996), The shape of Mars
+    and the topographic signature of the hemispheric dichotomy. Data columns
+    are "longitude," "latitude", and "radius (meters)."
+
+    The data are downloaded to a cache directory (usually ``~/.gmt/cache``) the
+    first time you invoke this function. Afterwards, it will load the data from
+    the cache. So you'll need an internet connection the first time around.
+
+    Returns
+    -------
+    data : pandas.DataFrame
+        The data table with columns "longitude", "latitude", and "radius(m)".
+    """
+    fname = which("@mars370d.txt", download="c")
+    data = pd.read_csv(
+        fname, sep="\t", header=None, names=["longitude", "latitude", "radius(m)"]
+    )
+    return data
diff --git a/pygmt/tests/test_datasets_samples.py b/pygmt/tests/test_datasets_samples.py
index ef598215822..1a52a3da2f3 100644
--- a/pygmt/tests/test_datasets_samples.py
+++ b/pygmt/tests/test_datasets_samples.py
@@ -6,6 +6,7 @@
     load_fractures_compilation,
     load_hotspots,
     load_japan_quakes,
+    load_mars_shape,
     load_ocean_ridge_points,
     load_sample_bathymetry,
     load_usgs_quakes,
@@ -76,6 +77,21 @@ def test_fractures_compilation():
     assert summary.loc["max", "azimuth"] == 360.0
 
 
+def test_mars_shape():
+    """
+    Check that the @mars370d.txt dataset loads without errors.
+    """
+    data = load_mars_shape()
+    assert data.shape == (370, 3)
+    summary = data.describe()
+    assert summary.loc["min", "longitude"] == 0.008
+    assert summary.loc["max", "longitude"] == 359.983
+    assert summary.loc["min", "latitude"] == -79.715
+    assert summary.loc["max", "latitude"] == 85.887
+    assert summary.loc["min", "radius(m)"] == -6930
+    assert summary.loc["max", "radius(m)"] == 15001
+
+
 def test_hotspots():
     """
     Check that the @hotspots.txt dataset loads without errors.