The present study focuses on the characterisation of the non-Gaussianity produced by radio point-sources and, for the first time, by infrared sources in the frequency range of the CMB from 30 to 350 GHz. We propose a simple prescription to infer the bispectrum from the power spectrum for point sources considering different independent populations of sources, with or without clustering. We test the accuracy of our prediction using publicly available all-sky simulations of radio and IR sources and find very good agreement. We further characterise the configuration dependence and the frequency behaviour of the IR and radio bispectra. We show that the IR bispectrum peaks for squeezed triangles and that the clustering of IR sources enhances the bispectrum values by several orders of magnitude at scales l ~ 100. At 150 GHz the bispectrum of IR sources starts to dominate over that of radio sources on large angular scales, and it dominates over the whole multipole range at 350 GHz. Finally, we compute the bias on fNL induced by radio and IR sources. We show that the positive bias induced by radio sources is significantly reduced by masking the sources. We also show, for the first time, that the form of the IR bispectrum mimics a primordial local bispectrum. The IR sources produce a negative bias which becomes important for Planck-like resolution and at high frequencies (Delta fNL ~ -6 at 277 GHz and Delta fNL ~ -60-70 at 350 GHz). Most of the signal being due to the clustering of faint IR sources, the bias Delta f_NL^IR is not reduced by masking sources above a flux limit and may, in some cases, even be increased due to the reduction of the shot-noise term.