That’s because I missed a zero. That works!
with pm.Model(coords={"time": df.index.values, "series": [0, 1]}) as m:
rho = pm.Normal("rho")
init = pm.Normal.dist()
sigma = pm.HalfNormal('sigma')
ar1 = pm.AR("ar1",
rho=rho,
sigma=sigma,
init_dist=init,
observed=df[['res1', 'res2']].to_numpy().T,
dims=("series", "time"),
shape=(2, 10000)
)