The connection between pharmacokinetic models and system theory has been established for a long time. In this approach, the drug concentration is seen as the output of a system whose input is the drug administered at different times. In this article we further explore this connection. We show that system theory can be used to easily accommodate any therapeutic regime, no matter its complexity, allowing the identification of the pharmacokinetic parameters by means of a non-linear regression analysis. We illustrate how to exploit the properties of linear systems to identify non-linearities in the pharmacokinetic data. We also explore the use of bootstrapping as a way to compare populations of pharmacokinetic parameters and how to handle the common situation of using multiple hypothesis tests as a way to distinguish two different populations. Finally, we demonstrate how the bootstrap values can be used to estimate the distribution of derived parameters, as can be the allometric scale factors.