How can q(z,l,u,c|x) be factorized to q(z|x)q(l|x)q(c|z)q(u|c,z)? How can q(z,l,u,c|x) be factorized to q(z|x)q(l|x)q(c|z)q(u|c,z)?
I don't know what assumptions are behind this equation, which is I want to know. I think normal chain rule of probability is in the process, but only chain rule can't derive the formula.
What I know is that z, l, u follows gaussian distribution, and c follows multinomial distribution.
It's from the paper Xu C, Lopez R, Mehlman E, Regier J, Jordan MI, Yosef N. Probabilistic harmonization and annotation of single-cell transcriptomics data with deep generative models. Mol Syst Biol. 2021 Jan;17(1):e9620. doi: 10.15252/msb.20209620. PMID: 33491336; PMCID: PMC7829634.
That formula is used for derive ELBO, since the model is similar to VAE.

In general $P(XY | A) = P(X | YA) P(Y | A)$.
If $X$ and $Y$ are conditionally independent given $A$ then $P(X | YA) = P(X | A) $ and $P(XY | A) = P(X | A) P(Y | A)$.
For your problem, if you assume the appropriate conditional independence structure for your variables, you can write that factorization.
I haven't looked at that paper, but something that's extremely common is to make some assumptions about what's conditionally independent of what to make the model more tractable and useful. In the best case, the assumptions are true (or in some sense close to true)! In other settings, the assumptions are even clearly false (e.g. Naive Bayes applied to language) but still produces useful enough results compared to alternatives that the assumptions are deemed worth making.
I often come back to the aphorism of George P. Box that "all models are wrong, but some are useful."
