The field of random vibrations of large-scale systems with millions of degrees-of-freedom (DOF) is of significant importance in many engineering disciplines. In this paper, we propose a method to calculate the time-dependent reliability of linear vibratory systems with random parameters excited by nonstationary Gaussian processes. The approach combines principles of random vibrations, the total probability theorem, and recent advances in time-dependent reliability using an integral equation involving the upcrossing and joint upcrossing rates. A space-filling design, such as optimal symmetric Latin hypercube (OSLH) sampling, is first used to sample the input parameter space. For each design point, the corresponding conditional time-dependent probability of failure is calculated efficiently using random vibrations principles to obtain the statistics of the output process and an efficient numerical estimation of the upcrossing and joint upcrossing rates. A time-dependent metamodel is then created between the input parameters and the output conditional probabilities allowing us to estimate the conditional probabilities for any set of input parameters. The total probability theorem is finally applied to calculate the time-dependent probability of failure. The proposed method is demonstrated using a vibratory beam example.