## Flat plane -- No interaction

x<-c(-10:10) # vector with numbers -10 to 10 in steps of 1

y<-c(-10:10) # vector with numbers -10 to 10 in steps of 1

f<-function(x,y){5+x+2*y} # the function I want to plot

# creates a matrix with values of f for each combination of values from x and y
z<-outer(x,y,f)

# the plot
persp(x,y,z, col="lightblue",theta=30, phi=20)


## Curved plane -- Including interaction term "x*y"

f<-function(x,y){5+x+2*y+x*y}
z<-outer(x,y,f)

persp(x,y,z, col="lightblue",theta=30, phi=20, 
             sub=expression(f(x,y)== 5+ x+2y))