Recently geographic information systems (GIS) have been widely used to analyze and preserve large amounts of data, text, and image in various fields related to physical geography, humanities and social sciences. Hotspot analysis is commonly applied in the spatial analysis of large amounts of data; however, the associated result is greatly determined by the sample method. Therefore, this study aims to explore the effect of grid size on the spatial distribution of temples in Taiwan through the hotspot analysis. The data source of this study was temple sites in the Government open data. The distribution of the top nine deities in Taiwan was obtained for related analysis. The hotspot distribution map with samples of different grid sizes was conducted by QGIS software. The geospatial characteristics of temples in Taiwan were analyzed through thematic maps such as point location maps, choropleth map, and the nearest neighbor analysis. The results show that the Nearest Neighbor Index (NNI) and Z values of the temples of top nine deities in the nearest neighbor analysis were all less than 1, which indicated an intense clustering. The Z values of Tudi Gong, Wangye had the mostly significant clustering with the range of -40 to -55, followed by those of Ma Zu, Sakyamuni, Guanyin, and Xuantian Shandi with the range of -30 to -40, and those of Holy Emperor Guan, Baosheng Emperor, and Zhongtan Yuanshuai with the range of -15 to -30. The ""expected average distance"" obtained from the nearest neighbor analysis was the most suitable radius for the performance of the heatmap. As the absolute value of Z is larger, the clustering distribution is more concentrated and the selection of grid size is more flexible. For example, the expected average distance of the Tudi Gong was about 5.8 km and the grid size with the radius of 7 km performed best in the hostspot distribution; however, the grid size with the radius of 5 km or 6 km also had similar distribution. In the distribution of the heatmap, the Tudi Gong, Holy Emperor Guan, and Sakyamuni are mainly distributed in the northern region, while Zhongtan Yuanshuai, Wangye, Xuantian shandi, Baosheng Emperor, Ma Zu, and Guanyin are distributed in the southern region. Among them, Wangye and Baosheng Emperor are particularly concentrated in the western coastal areas, while Xuantian shandi and Ma Zu are distributed in the western hilly areas. The results of kernel density estimation and the Getis-OrdGi algorithm were slightly different in the distribution of hotspot areas. It suggests that the nearest neighbor analysis should be firstly conducted to obtain the expected average distance, then the grid size with optimal radius could be determined based on the Z values for further study.