Hi Junpeng, sorry for the late response. I have some difficulties sharing my code and hope it would work for me to describe what I did.
The model itself is following document: https://docs.pymc.io/notebooks/LKJ.html (same hyper-parameters) assuming data is from multivariate normal distribution with dims = (10000,3), where my mu = [1,10,-5], cov=[[1,0.1,0.2],[0.1,1,0.9],[0.2,0.9,1]]. Data was then masked with np.nan each column with missing rate = 0.3. I used ADVI as inference. It only took minutes to get the trace and interesting things started happening from here - not only well recovered cholesky matrix but give good predictions on missing