Normalization to reference genes is the most common method to avoid bias in real-time quantitative PCR (qPCR), which has been widely used for quantification of gene expression. Despite several studies on gene expression, Lilium, and particularly L. regale, has not been fully investigated regarding the evaluation of reference genes suitable for normalization. In this study, nine putative reference genes, namely 18S rRNA, ACT, BHLH, CLA, CYP, EF1, GAPDH, SAND and TIP41, were analyzed for accurate quantitative PCR normalization at different developmental stages and under different stress conditions, including biotic (Botrytis elliptica), drought, salinity, cold and heat stress. All these genes showed a wide variation in their Cq (quantification Cycle) values, and their stabilities were calculated by geNorm, NormFinder and BestKeeper. In a combination of the results from the three algorithms, BHLH was superior to the other candidates when all the experimental treatments were analyzed together; CLA and EF1 were also recommended by two of the three algorithms. As for specific conditions, EF1 under various developmental stages, SAND under biotic stress, CYP/GAPDH under drought stress, and TIP41 under salinity stress were generally considered suitable. All the algorithms agreed on the stability of SAND and GAPDH under cold stress, while only CYP was selected under heat stress by all of them. Additionally, the selection of optimal reference genes under biotic stress was further verified by analyzing the expression level of LrLOX in leaves inoculated with B. elliptica. Our study would be beneficial for future studies on gene expression and molecular breeding of Lilium.