The analytical solution for seismic wave propagation associated with a point force in a fluid‐filled porous medium is developed. The point force solution is applied to solve the boundary value problem of seismic wave propagation in a stratified poroelastic medium. The coupled Biot vector wave equations are expressed in cylindrical coordinates and expanded in Fourier series with respect to azimuth. The resulting equations are transformed to the wave‐number domain using the Hankel transform method. Following this analysis, one set of three coupled partial differential equations associated with fast compressional, slow compressional, and vertically polarized shear waves is derived. The unknowns are the radial and vertical displacements associated with the solid frame motion as well as the fluid pore pressure. A separate partial differential equation associated with waves whose particle motion is polarized in horizontal planes (SHwaves) is derived as well. The general solution of the three coupled differential equations is obtained by the Kupradze method.This solution of the Biot’s motion equations in the wave‐number domain leads directly to closed form expressions for the vector wave displacement and the pressure produced by a point force in a poroelastic unbounded medium. In order to develop the solution of a point force in the presence of a stratified porous medium, the displacement‐stress matrix, the pressure, and the vertical component of the displacement of the fluid in its relative motion versus the solid for different regions are expressed in terms of upgoing and downgoing waves and unknown wave coefficients. The wave coefficients are determined by applying boundary conditions of continuity of displacements, pressure and stresses across each layer interface, and the radiation conditions at infinity. To determine the unknown wave coefficients, a method that consists in expressing the kernels of the Hankel transform integrals in terms of factorization of upgoing and downgoing wave amplitudes in each layer is used. These factorizations are based on the generalized reflection and transmission coefficient matrices, which are formed recursively, from one layer boundary to the next, including all of the reflection/conversion/transmission properties of the layered medium.Their factorization method allows the field within each layer above or below the source to be determined once the field in the medium containing the source is known. The final equations provide a complete description of the field throughout the layered medium. The particular form of the equations makes possible the simultaneous evaluation of the response at a number of detector locations for a number of different source positions in a borehole for interwell seismic applications. Numerical model results demonstrate the validity of this theoretical development for predicting spectral responses associated with porosity and permeability effects. The seismic pressure response of a thin gas‐saturated porous layer was analyzed. The results inferred that the gas‐saturated porous layer strongly attenuates the waveforms observed by detectors within the layer. Alternatively, large multiple reflections and converted waves from the layer are observed by detectors in the water‐saturated porous formation.