Otdfbt

Look at this image:

This is what I get from this:

begintikzpicture
 draw[style=help lines] (-5,-5) grid (5,5);
 draw (-4,0)--(4,0);
 draw (0,-4)--(0,4);
 foreach y in -4,-3,...,4 
 draw (0 - 0.1,y) -- (0+0.1,y);
 draw (y,0 - 0.1) -- (y,0+0.1);
 

 %Nodes:
 node (a0) at (-4,-4) ;
 draw[fill] (a0) circle [radius=1.5pt];
 node (a1) at (-2,4) ;
 draw[fill] (a1) circle [radius=1.5pt];
 node (a2) at (2,-2) ;
 draw[fill] (a2) circle [radius=1.5pt];
 node (a3) at (4,2) ;
 draw[fill] (a3) circle [radius=1.5pt]; 

 draw (-4,-4) to (-2,4) to (2,-2) to (4,2); % to (a2) to (a3);

endtikzpicture

I'm trying to to get a line between them (the dots) that will be like a function (not a straight line - a curve like a polynomial).

Is this possible?

Thank you!

You can use plot [smooth] coordinates (which is not a single polynom but a spline):

documentclass[tikz]standalone

begindocument
begintikzpicture
 draw[style=help lines] (-5,-5) grid (5,5);
 draw (-4,0)--(4,0);
 draw (0,-4)--(0,4);
 foreach y in -4,-3,...,4 
 draw (0 - 0.1,y) -- (0+0.1,y);
 draw (y,0 - 0.1) -- (y,0+0.1);
 

 %Nodes:
 node (a0) at (-4,-4) ;
 draw[fill] (a0) circle [radius=1.5pt];
 node (a1) at (-2,4) ;
 draw[fill] (a1) circle [radius=1.5pt];
 node (a2) at (2,-2) ;
 draw[fill] (a2) circle [radius=1.5pt];
 node (a3) at (4,2) ;
 draw[fill] (a3) circle [radius=1.5pt]; 

 draw plot [smooth] coordinates (-4,-4) (-2,4) (2,-2) (4,2); % to (a2) to (a3);

endtikzpicture
enddocument

Solution which forces the middle points to have a horizontal tangent:

documentclass[tikz,border=3.14]standalone

begindocument
begintikzpicture
 draw[style=help lines] (-5,-5) grid (5,5);
 draw (-4,0)--(4,0);
 draw (0,-4)--(0,4);
 foreach y in -4,-3,...,4 
 draw (0 - 0.1,y) -- (0+0.1,y);
 draw (y,0 - 0.1) -- (y,0+0.1);
 

 %Nodes:
 node (a0) at (-4,-4) ;
 draw[fill] (a0) circle [radius=1.5pt];
 node (a1) at (-2,4) ;
 draw[fill] (a1) circle [radius=1.5pt];
 node (a2) at (2,-2) ;
 draw[fill] (a2) circle [radius=1.5pt];
 node (a3) at (4,2) ;
 draw[fill] (a3) circle [radius=1.5pt]; 

 draw (-4,-4) to[out=90,in=180] (-2,4) to[out=0,in=180] (2,-2) to[out=0,in=-95] (4,2); % to (a2) to (a3);

endtikzpicture
enddocument

I don't know how to compute this in LaTeX easily, so I fitted a plot using Python's numpy.polyfit and used the result to plot the fit in TikZ:

documentclass[tikz,border=3.14]standalone

%% polynomial coefficients found with Python (numpy.polyfit)
%% $f(x) = 0.1875 x^3 - 1/6 x^2 - 2.25 x^1 + 10/6 x^0$

begindocument
begintikzpicture
 draw[style=help lines] (-5,-5) grid (5,5);
 draw (-4,0)--(4,0);
 draw (0,-4)--(0,4);
 foreach y in -4,-3,...,4 
 draw (0 - 0.1,y) -- (0+0.1,y);
 draw (y,0 - 0.1) -- (y,0+0.1);
 

 %Nodes:
 node (a0) at (-4,-4) ;
 draw[fill] (a0) circle [radius=1.5pt];
 node (a1) at (-2,4) ;
 draw[fill] (a1) circle [radius=1.5pt];
 node (a2) at (2,-2) ;
 draw[fill] (a2) circle [radius=1.5pt];
 node (a3) at (4,2) ;
 draw[fill] (a3) circle [radius=1.5pt]; 

 draw plot[domain=-4:4,samples=100] (x, .1875*x*x*x - x*x/6 - 2.25*x + 10/6);

endtikzpicture
enddocument

Just for your information. You can calculate and plot the interpolation polynomial with Python and the two libraries Matplotlib and NumPy:

import numpy as np
import matplotlib.pyplot as plt

x = (-4, -2, 2, 4)
y = (-4, 4, -2, 2)
p = np.polyfit(x,y,3)
t = np.linspace(min(x),max(x),num=100)
f = np.polyval(p,t)
plt.plot(t,f)

Matplotlib supports export to TikZ code (actually it exports to PGF) and to save the plots directly as PDF created with TikZ and LaTeX (see for example https://tex.stackexchange.com/a/426071/117050 and https://tex.stackexchange.com/a/391078/117050 for some code that might get you started).

We can use draw controls - the red curve, in comparison with the blue curve draw plot[smooth] coordinates. (if you want, you can control so that the red and blue curves are identical)

documentclass[tikz,border=5mm]standalone
begindocument
begintikzpicture
draw[gray!30] (-5,-5) grid (5,5);
draw (-5,0)--(5,0) (0,-5)--(0,5);
foreach i in -5,...,5
draw 
(0,i)--+(1mm,0)--+(-1mm,0)
(i,0)--+(0,1mm)--+(0,-1mm);

draw[blue] plot[smooth] coordinates 
(-4,-4) (-2,4) (2,-2) (4,2);

draw[red] 
(-4,-4)..controls +(80:1) and +(180:1)..
(-2,4)..controls +(0:1) and +(180:1)..
(2,-2)..controls +(0:1) and +(-100:1)..
(4,2);

foreach p in (-4,-4),(-2,4),(2,-2),(4,2)
fill p circle(2pt);
endtikzpicture
enddocument 

With some calculations, I found formula of the function is
-(1/72)*x^4+3/16*(x^3)+(1/9)*x^2-9/4*x+7/9
I use pgfplots to draw

documentclass[tikz]standalone
usepackagepgfplots
pgfplotssetcompat=1.16
usepackagefouriernc
begindocument
 begintikzpicture[
 declare function=
 f(x)=-(1/72)*x^4+3/16*(x^3)+(1/9)*x^2-9/4*x+7/9;
 
 ]
 beginaxis[axis equal,
 width=12 cm, 
 grid=major, 
 axis x line=middle, axis y line=middle,
 axis line style = very thick,
 grid style=gray!30,
 ymin=-5, ymax=5, yticklabels=, ylabel=$y$,
 xmin=-5, xmax=5, xticklabels=, xlabel=$x$,
 samples=500,
 ]
 addplot[blue, very thick,domain=-5:5, smooth]f(x);
 node[below] at (-2, 0) $-2$;
 node[above ] at (-4, 0) $-4$;
 node[below ] at (4, 0) $4$;
 node[right] at (0,-4) $-4$;
 node[left ] at (0,2) $2$;
 node[ right ] at (0,4) $4$;
 node[below right] at (0, 0) $O$;
 node[above ] at ( 2,0) $2$;
 node[left ] at (0, -2) $-2$;
 addplot [mark=*,only marks,samples at=-4,-2,2,4] f(x);
 ;
 draw[dashed, thick] (-4,0) -- (-4,-4) -- (0,-4);
 draw[dashed, thick] (-2,0) -- (-2,4) -- (0,4);
 draw[dashed, thick] (2,0) -- (2,-2) -- (0,-2);
 draw[dashed, thick] (4,0) -- (4,2) -- (0,2);
 endaxis
 endtikzpicture
enddocument 

Results from Maple.

With marmot's help , I reduce my code

documentclass[tikz]standalone 
usepackagepgfplots 
pgfplotssetcompat=1.16 
usepackagefouriernc 
begindocument 
begintikzpicture[ 
declare function= 
f(x)=-(1/72)*x^4+3/16*(x^3)+(1/9)*x^2-9/4*x+7/9; 
 
] 
beginaxis[axis equal, 
width=12 cm, 
grid=major, 
axis x line=middle, axis y line=middle, 
axis line style = very thick, 
grid style=gray!30, 
ymin=-5, ymax=5, yticklabels=, ylabel=$y$, 
xmin=-5, xmax=5, xticklabels=, xlabel=$x$, 
samples=500, 
] 
addplot[blue, very thick,domain=-5:5, smooth]f(x); 

addplot [mark=*,only marks,samples at=-4,-2,2,4] f(x); 
; 

pgfplotsinvokeforeach-4,-2,2,4draw[dashed] (#1,0) 

foreach X/Y in -4/right,-2/left,2/left,4/right 
edeftempnoexpandnode[Y] at (0,X) $X$; 
temp 


foreach X/Y in -4/above,-2/below,2/above,4/below 
edeftempnoexpandnode[Y] at (X,0) $X$; 
 temp 
% 
endaxis 
endtikzpicture 
enddocument