The B mode of the polarization of the cosmic microwave background, which is a unique estimator of the presence of cosmological gravitational waves in the early universe, is the target of upcoming high sensitivity CMB experiments. In this thesis, new methods are developed for statistical analysis of the CMB polarization. The E and B modes are defined in the space of the Stokes parameters, and the mathematical properties of the E-B decomposition are investigated. These E and B families are applied to the analysis of the morphology and frequency dependence of the thermal dust polarized foregrounds, revealing spectral differences between the E and B modes. This is followed by the presentation of a new method for correction of E-B leakage and the estimation of the B mode map using partial sky data. The method can improve the accuracy of the B mode power spectrum estimation by multiple orders of magnitude. The statistics of the polarization angle and the non-polarized points are investigated. Among other results, the E and B modes are found to have distinct signatures in the polarization angle distribution function, and the number density of non-polarized points is shown to correlate with the tensor-to-scalar ratio under conditions of gravitational lensing. These statistical methods are applied to the polarization data from the latest Planck release. Lastly, the large-scale anomalies of the CMB are examined using new map-space methods, which reveal connections in the parity asymmetry and the alignments of the low multipoles.