A standard 1-compartment TK model has 2 parameters, so you cannot change the model in BYOM to make it 3. You have a choice which 2 parameters to use (ke and ku, or ke and BCF, or even ku and BCF). The third variable will then be a model output. So, technically, you would need to propagate the uncertainty in 2 parameters to a model output. Parts of BYOM do that for specific output such as the LCx,t. However, there is no facility for BCF, nor any general-purpose facility (yet).
It is possible to add code to read the sample from the saved MAT file and code an ad-hoc error propagation. However, that requires quite a bit of BYOM-savvy, and can easily go wrong. The way you did it (make two model versions) certainly works, but it is a bit of double work.
There is a functionality in BYOM for using a zero-variate data point: a data point without a time vector, which is a translation of the model parameters. However, at this point, the model prediction for this point is not plotted with a CI. I explicitly excluded that in the code, and I am not sure why

. I will put it on my to-do list to look into this. As soon as I have a bit of time I can see if this functionality in BYOM can be epxanded.