Function that is a good fit to the plot made with SmoothHistogram

SmoothHistogram uses kernel density estimation. So the correct and direct way to obtain the function is to construct a SmoothKernelDistribution, then take its PDF.

distr = SmoothKernelDistribution[data]

Plot[PDF[distr, x], {x, 0, 6}]

This is one approach yielding an interpolating function:

  data = {4, 4, 4, 4, 1, 1, 4, 1, 5, 5, 5, 5, 5};
  f = Interpolation[
     DeleteDuplicates[(List @@ 
          First@First@
             Cases[ SmoothHistogram[data] , _Line, 
                      Infinity]), #1[[1]] == #2[[1]] &]]

  Plot[f[x], {x, 0, 6.5}]