Volcanic activity often results in the formation of new volcanic vents. These new vents can create hazards in unexpected areas. Therefore, the probability of new vent formation should be assessed as part of volcanic hazard assessments. This paper describes our use of kernel density estimation (KDE) as a way to estimate the spatial density of future volcanic vents. The bivariate Gaussian kernel function is described step-by-step using pseudocode. Our computer code, written in PERL, is used to calculate the spatial density of existing vents and then create a contour map using GMT (Generic Mapping Tools). Application of this method and code relies on several assumptions about the definition of volcanic events, independence of events, the type of kernel function used, and the selection of kernel bandwidth. Three examples using the code are provided: (1) for volcanic vents located west of the city of Managua (Nicaragua), (2) for volcanic vents distributed within the Arsia Mons caldera (Mars), weighted by volume, and (3) for vents of the Lassen volcanic system (northern California), sub-divided by geochemistry.