Bolt connection is one of the main fixing methods of cylindrical shell structures. A typical bolted connection model is considered as a tuned system. However, in the actual working conditions, due to the manufacturing error, installation error and uneven materials of bolts, there are always random errors between different bolts. To investigate the influence of non-uniform parameters of bolt joint, including the stiffness and the distribution position, on frequency complexity characteristics of cylindrical shell through a statistical method is the main aim of this paper. The bolted joints considered here were simplified as a series of springs with random features. The vibration equation of the bolted joined cylindrical shell was derived based on Sanders’ thin shell theory. The Monte Carlo simulation and statistical theory were applied to the statistical analysis of mode characteristics of the system. First, the frequency and mode shape of the tuned system were investigated and compared with FEM. Then, the effect of the random distribution and the random constraint stiffness of the bolts on the frequency and mode shape were studied. And the statistical analysis on the natural frequencies was evaluated for different mistuned levels. And some special cases were presented to help understand the effect of random mistuning. This research introduces random theory into the modeling of bolted joints and proposes a reference result to interpret the complexity of the modal characteristics of cylindrical shells with non-uniform parameters of bolt joints.