[GIS] How to superimpose a geopandas dataframe on a raster plot

Question

I have a geopandas dataframe countries with country polygons, and a raster dataset raster (read with rasterio). I've already aligned the two as follows:

countries = countries.to_crs(raster.crs.data)

I can see they're aligned when I plot each separately. However, I have been unable to display the country borders overlaid on top of the raster. The following just shows the raster image, no country borders:

import matplotlib.pyplot as plt    
from rasterio import plot as rioplot
plt.figure(figsize=(15,15))
base = rioplot.show(raster)
countries.plot(ax=base, color='red');

Answer 1

For me it works if I pass the matplotlib ax object explicitly to rasterio.plot.show:

fig, ax = plt.subplots(figsize=(15, 15))
rasterio.plot.show(raster, ax=ax)
countries.plot(ax=ax, facecolor='none', edgecolor='r');

Full example with raster data from naturalearth (link) (both datasets have the same CRS, so I don't do any reprojection):

import geopandas
import matplotlib.pyplot as plt
import rasterio
import rasterio.plot

countries = geopandas.read_file(geopandas.datasets.get_path('naturalearth_lowres'))
raster = rasterio.open("NE1_50M_SR_W/NE1_50M_SR_W.tif")

fig, ax = plt.subplots(figsize=(15, 15))
rasterio.plot.show(raster, ax=ax)
countries.plot(ax=ax, facecolor='none', edgecolor='red')

gives:

For me it also does not work without passing ax. Not sure whether this is a bug in rasterio or geopandas.