Background: Spatial and temporal resolution of brain network activity can be improved by combining different modalities. Functional Magnetic Resonance Imaging (fMRI) provides full brain coverage with limited temporal resolution, while electroencephalography (EEG), estimates cortical activity with high temporal resolution. Combining them may provide improved network characterization. New Method: We examined relationships between EEG spatiospectral pattern timecourses and concurrent fMRI BOLD signals using canonical hemodynamic response function (HRF) with 1st and 2nd temporal derivatives in voxel-wise general linear models (GLM). HRF shapes were derived from EEG-fMRI time courses during "resting-state", visual oddball and semantic decision paradigms. Results: The resulting GLM F-maps self-organized into several different large-scale brain networks (LSBNs) often with different timing between EEG and fMRI revealed through differences in GLM-derived HRF shapes (e.g., with a lower time to peak than the canonical HRF). We demonstrate that some EEG spatiospectral patterns (related to concurrent fMRI) are weakly task-modulated. Comparison with existing method(s): Previously, we demonstrated 14 independent EEG spatiospectral patterns within this EEG dataset, stable across the resting-state, visual oddball and semantic decision paradigms. Here, we demonstrate that their time courses are significantly correlated with fMRI dynamics organized into LSBN structures. EEG-fMRI derived HRF peak appears earlier than the canonical HRF peak, which suggests limitations when assuming a canonical HRF shape in EEG-fMRI. Conclusions: This is the first study examining EEG-fMRI relationships among independent EEG spatiospectral patterns over different paradigms. The findings highlight the importance of considering different HRF shapes when spatiotemporally characterizing brain networks using EEG and fMRI.