Molecular communication (MC) is recently featured as a novel communication tool to connect individual biological nanorobots. It is expected that a large number of nanorobots can form large multi-agent MC systems through MC to accomplish complex and large-scale tasks that cannot be achieved by a single nanorobot. However, most previous models for MC systems assume a unidirectional diffusion communication channel and cannot capture the feedback between each nanorobot, which is important for multi-agent MC systems. In this paper, we introduce a system theoretic model for large-scale multi-agent MC systems using transfer functions, and then propose a method to analyze the stability for multi-agent MC systems. The proposed method decomposes the multi-agent MC system into multiple single-input and single-output (SISO) systems, which facilitates the application of simple analysis technique for SISO systems to the large-scale multi-agent MC system. Finally, we demonstrate the proposed method by analyzing the stability of a specific large-scale multi-agent MC system and clarify a parameter region to synchronize the states of nanorobots, which is important to make cooperative behaviors at a population level.