## 2 Methodology: Null Space on Scaled Independent Variable  
The goal of the method presented in this section is to fit discrete data measurements \((s,\ y_{s})\) with a real valued rational function \(f\) where \(s\) is the explanatory variable and \(y_{s}\) is the response variable. We consider the fitting function \(f\) to be in the form  
\[f(z)=\frac{a_{0}+a_{1}z+a_{2}z^{2}+\cdots a_{N}z^{N}}{b_{0}+b_{1}z+b_{2}z^{2}+ \cdots b_{N}z^{N}} \tag{1}\]  
for real value \(z\). When the responses \(y_{s}\) are complex numbers, as in many engineering systems, the fitting technique can be broken down into fitting separately the real and imaginary parts of \(y_{s}\). Before we proceed with the derivations, we normalize the explanatory variable data \(s\) in the infinity norm sense. So, all values \(s\) will be replaced with \(\frac{s}{s_{max}}\), where \(s_{max}\) is the maximum value among all \(\left\lvert s\right\rvert\). This scaling process allow all powers of \(s\) to be contained in the interval \([0,1]\) for positive \(s\). For each value \(s\), Equation (1) becomes:  
\[y_{s}(b_{0}+b_{1}s+b_{2}s^{2}+\cdots b_{N}s^{N})-(a_{0}+a_{1}s+a_{2}s^{2}+ \cdots a_{N}s^{N})=0\]which gives  
\[b_{0}y_{s}+b_{1}sy_{s}+b_{2}s^{2}y_{s}+\cdots b_{N}s^{N}y_{s}-a_{0}-a_{1}s-a_{2}s^ {2}-\cdots-a_{N}s^{N}=0\]  
The unknowns are the coefficients in the vector \((b_{0},b_{1},b_{2},\ldots,b_{N},a_{0}\,,a_{1}\,,a_{2}\ldots a_{N})\). Once we solve for the coefficients, a fitting function is defined. Suppose we have m data (\(s\), \(y_{s}\)). Each pair of measurement gives an equation, and therefore we obtain a homogeneous system of m equations represented in matrix form:  
\[\mathbf{Ax}=0 \tag{2}\]  
where \(\mathbf{x}=(b_{0},b_{1},b_{2},\ldots,b_{N},a_{0}\,,a_{1}\,,a_{2}\ldots a_{N})\)' and A is the coefficient matrix with size m by \(2N+2\). Each couple (\(s\), \(y_{s}\)) defines a row in the matrix A,  
\[(y_{s}\;\;y_{s}s\;\;y_{s}s^{2}\ldots\;y_{s}s^{N}\;-1\;-s\;\;-s^{2}\ldots-s^{N}\;)\]  
The benefit of scaling by \(1/s_{max}\), where all powers of s are contained in the interval \([0,1]\), is evident. This avoids a badly scaled coefficient matrix A, especially for higher \(s\) and \(N\). The set of vector solutions \(\mathbf{x}=(\mathbf{x}_{1},\mathbf{x}_{2},\ldots,\mathbf{x}_{r})\) form a basis of the null space of A, null(A). With the proper choice of the degree \(N\), each solution vector \(\mathbf{x}_{i}\,,\;1\leq i\leq r\) gives a rational fitting function \(f_{i}\). The best rational function \(f\) is picked by calculating the minimum of the errors \(\left.||\mathbf{f}_{i}-\mathbf{y}_{s}||_{\infty}\right.\), where \(\mathbf{f}_{i}\) is the evaluation vector with inputs \(s\) using Equation (1), and \(\mathbf{y}_{s}\) is the vector of the data responses \(y_{s}\). In other words, the index n, that gives the best fitting function \(f_{n}\) is defined by  
\[||\mathbf{f}_{n}-\mathbf{y}_{s}||_{\infty}=\min_{1\leq i\leq r}||\mathbf{f}_{ i}-\mathbf{y}_{s}||_{\infty}\]  
The vector \(\mathbf{f}_{n}\) contains the coefficients of the model function \(f_{n}\). This process is our preferred error minimization method between the response observations \(\mathbf{y}_{s}\) and the function coefficients in \(\mathbf{f}_{i}\). The users can choose an alternative minimization norm such as the 2-norm for instance. Next, we horizontally stretch the fitting function by a factor of \(s_{max}\) to be defined on the original interval containing the raw unscaled variable \(s\), i.e., the fitting function is giving by  
\[f(z)=f_{n}(s_{max}z) \tag{3}\]  
To solve Equation (2) we used MATLAB built-in function null(A), which is based on the singular value decomposition (SVD) of matrix A. When A is very large the economy-sized MATLAB svd(A,\"econ\") can be used. This method also extends for negative values of \(s\), in which case scaled values \(\frac{s}{s_{max}}\) will be contained in the interval \([-1,1]\), and so are all powers of s. The existence of infinitely many solutions is guaranteed when the system (2) is underdetermined \((m<2N+2)\).