simulate and estimate ARMA model

I'm trying to simulate an ARMA model and later simulate the AR and MA polynomial coefficients from the simulated time series.

I found the following ARIMA function in MatLab that I hoped help me achieve my goal. Looking at the provided examples, I thought this would work:

AR_coef = {0.5, -0.1};
MA_coef = {0.8, 0.5, 0.3};

% simulate time series
Mdl_in = arima('AR', AR_coef, ...
    'MA', MA_coef, ...
    'Constant', 0, ...
    'Variance', 1);
y = simulate(Mdl_in,1,'NumPaths',10000);

% estimate coefficients 
Mdl_out = arima(numel(AR_coef), 0, numel(MA_coef));
T = numel(y);
idxpre = 1:Mdl_out.P;
idxest = (Mdl_out.P + 1):T;
EstMdl = estimate(Mdl_out,y(idxest)','Y0',y(idxpre)');

My hope was that EstMdl returns the AR_coef and MA_coef but it does not. Based on the random number generation, the estimated model parameter varies a lot. Can someone help me generate an example where I simulate an ARMA model and later estimate its coefficients?

I'm also fine with using Python if it is easier.