A point collocation nonintrusive polynomial chaos approach is presented for nonlinear aeroelastic uncertainty quantification analysis. A unique aspect of the approach is the use of a nonlinear frequency-domain harmonic balance aeroelastic solution methodology, which allows one to directly solve for the data points that comprise an aeroelastic Limit Cycle Oscillation (LCO) unsteady first harmonic pitching motion response amplitude versus reduced velocity curve. Other researchers have typically framed their aeroelastic uncertainty quantification analysis methods around time-domain aeroelastic solvers, where an indirect measure of aeroelastic stability, such as aeroelastic damping, is the dependent variable for which one performs an uncertainty quantification analysis. The presented method allows one to directly determine uncertainty quantification for an aeroelastic limit cycle oscillation response curve, including possible unstable branches of such curves, which would be challenging with a time domain approach. The technique is demonstrated for uncertainty in freestream Mach number and airfoil structural static unbalance for a two-dimensional model of the NASA Benchmark Models Program B0012 NACA 0012 airfoil section aeroelastic wing configuration.