The spatio‐temporal variability of aerosols over Indian subcontinent is mainly due to transported dust from adjacent deserts (Thar and Middle East deserts), local emission due to anthropogenic activities and prevailing meteorological conditions. On large scale, the quantification of transported and locally emitted aerosol from these regions is complicated. Here, we used empirical orthogonal function (EOF) analysis to identify these regions and their variability by using 33 years of Total Ozone Monitoring Spectrometer (TOMS) satellite data. The maximum variability in aerosol is explained by first two EOF modes (70.31 and 20.57%) over Indian subcontinent. The major aerosol, i.e. transported dust from adjacent deserts confined to NW India and Pakistan, is observed in first leading mode, whereas biomass burning, industrial and dense populated region of southeast and eastern region of the Indo‐Gangetic Plain (IGP) are revealed in the second dominant mode. The EOF analysis is carried out specifically for the pre‐monsoon and monsoon seasons over Indian subcontinent as maximum aerosol loading is observed during this period. The region of NW India, IGP, Pakistan and northern Arabian Sea explains maximum variability in both seasons. The first three leading modes and their relationship with different atmospheric and surface variables are carried out for pre‐monsoon and monsoon seasons. This study explains the potential role of aerosols on reduction in cloudiness, increased shortwave at the ground, land‐surface tropospheric warming and its feedback to other related processes. This study strongly suggests that there is a need for further appropriate observational as well as modelling study on the role of semi‐direct aerosol effect over Indian subcontinent.