Following the decomposition of the correlation function in Eq. ( ), one can introduce a series of auxiliary density operators (ADOs), denoted as ρ n (p t c ,x t c ,t). Each ADO is specified by a superindex n = (⋯,n α p ,⋯), where n α p is a non-negative integer and can be assigned as the occupation number of the pth effective dissipative boson in bath α. The equation of motion for the ADOs is given by