How to approach to fitting curve? I'd like to approximate fitting curve some kind of curves like below.
(1, 3.5), (2, 4.3), (3, 7.2), (4, 8) which is having 4 points.
 
and I heard that this solver is PINV() of matlab function.
But I don't know how to use.
Would you please let me know how to use and find a approximation fitting curve equation?

Problem statement
Fit the given data set with a sequence of polynomial fits:
$$
 y(x) = a_{0} + a_{1} x + \dots + a_{d} x^{d}
$$
Indicate where we can solve the linear system with the normal equations, and when we must rely exclusively upon the pseudoinverse.
Linear System
The linear system for the polynomial with highest order $d$ is
$$
\begin{align}
  \mathbf{A} a &= y \\
% A
\left[ \begin{array}{cccc}
  1 & x_{1} & \cdots & x_{1}^{d} \\
  1 & x_{2} & \cdots & x_{2}^{d} \\
  1 & x_{3} & \cdots & x_{3}^{d} \\
  1 & x_{4} & \cdots & x_{4}^{d} \\
\end{array} \right]
% a
\left[ \begin{array}{cccc}
  a_{0} \\
  a_{1} \\
  \vdots \\
  a_{d} \\
\end{array} \right]
%
&=
%
% y
\left[ \begin{array}{c}
  y_{1} \\
  y_{2} \\
  y_{3} \\
  y_{4} 
\end{array} \right]
%
\end{align}
$$
Least squares solution
The least squares solution is defined as
$$
 a_{LS} = \left\{
a \in \mathbb{R}^{d+1} \colon
\lVert
\mathbf{A} a - y
\rVert_{2}^{2}
\text{ is minimzed}
\right\}
$$
The residual error vector is
$$
 r = \mathbf{A} a - y.
$$
The total error, the quantity minimized, is $r^{2} = r\cdot r$.
Summary of results
The total error demonstrates typical behavior. Increasing the order of fit initially reduces the error. Then it plateaus, or may actually increase. 
For these data, the cubic fit provides  the best combination of total error and computational cost.

The amplitudes for each order are collected in the table below.
$$
\begin{array}{clllllll}
 k & a_{0} & a_{1} & a_{2} & a_{3} & a_{4} & a_{5} & a_{6} \\
 0 & 5.75 \\
 1 & 1.65 & \phantom{-}1.64 \\
 2 & 1.65 & \phantom{-}1.64 & \phantom{-}0  \\
 \color{blue}{3} & \color{blue}{9.} & \color{blue}{-10.05} & \phantom{-}\color{blue}{5.25} & \color{blue}{-0.7} \\
 4 & 4.04271 & \phantom{-}0.277692 & -1.97938 & \phantom{-}1.36554 &
   -0.206554 \\
 5 & 2.57594 & \phantom{-}1.38279 & -0.0545302 & -0.868031 &
   \phantom{-}0.545109 & -0.0812778 \\
 6 & 1.89099 & \phantom{-}1.41016 & \phantom{-}0.686033 & -0.172893 &
   -0.615848 & \phantom{-}0.350573 & -0.0490168 \\
\end{array}
$$
Lower order fits: $d\le 3$
For $d=0$, the amplitude is given by the mean of the data $a_{0} = \bar{y}$.
For $d=1,2,3$, the normal equations
$$
 \mathbf{A}^{T} \mathbf{A} a = \mathbf{A}^{T} y
$$
can be solved as
$$
  a_{LS} = \left( \mathbf{A}^{T} \mathbf{A} \right)^{-1} \mathbf{A}^{T} y
$$
Higher order fits: $d\ge 4$
For $d\ge4$, the product matrix $\mathbf{A}^{T} \mathbf{A}$ is rank deficient and cannpt be inverted. The solution requires the pseudoinverse. The least squares minimizers are the affine space
$$
 a_{LS} = \mathbf{A}^{+}y + \left( \mathbf{I}_{4} - \mathbf{A}^{+} \mathbf{A} \right) z, \quad  z\in \mathbb{R}^{4}
$$
Sequence of plots

