Design a site like this with WordPress.com
Get started

Computing a competitive equilibrium with multiple goods in R

Simple visuals of something remarkable.

In the previous blog post titled “Computing the demand for multiple goods in R.” I discussed how we can use embedded for loops along with the command solnl() from the NlcOptim package in r to derive multiple demand curves and make comparisons of cross price changes.

What we will do now is add a supply side to our model where we will consider two different firms which act as price takers in two distinct markets.

When modelling the firms we will employ a different optimization function that being optim() which is a base R function. The reason why we use this over solnl() is because we don’t need to “trick” the optimizer into doing what we want it to do by reframing the profit maximization problem as minimizing a negative profit function. In this case we just use a built in command within the function.

Lets see the code in which we generate our data.

library('NlcOptim')
library('ggplot2')
library('ggpubr')
library('ggthemes')

#set up model parameters
p<-seq(1,20,by=1)
m<-100

#storage for demands,supply and prices
x1d<-NULL 
x2d<-NULL
p1data<-NULL
p2data<-NULL
x1s<-NULL
x2s<-NULL

#Guess
guess<-c(1,1)

#Utility function
pref<-function(x){
  return(-x[1]^0.5*x[2]^0.5)
}

#constraint
con<-function(x){
  f=NULL
  f=rbind(f,p1*x[1]+p2*x[2]-m)
  return(list(ceq = NULL, c = f))
}

#profit functions
firm1<-function(x){
  return(p1*x-2*x^2-2*x)
}

firm2<-function(x){
  return(p2*x-x^2-x)
}




#Loop for Demand
for(p1 in p){
  for(p2 in p){
    umax<-solnl(guess,pref,con)
    d1<-umax$par[1,1]
    d2<-umax$par[2,1]
    x1d<-rbind(x1d,d1)
    x2d<-rbind(x2d,d2)
    p1data<-rbind(p1data,p1)
    p2data<-rbind(p2data,p2)
  }
}

#Loop for Supply
for(p1 in p){
  for(p2 in p){
    pmax1<-optim(guess[1],firm1,control=list(fnscale=-1),method='L-BFGS-B',lower=0)
    pmax2<-optim(guess[1],firm2,control=list(fnscale=-1),method='L-BFGS-B',lower=0)
    s1<-pmax1$par[1]
    s2<-pmax2$par[1]
    x2s<-rbind(x2s,s2)
    x1s<-rbind(x1s,s1)
  }
}

#create excess demand functions
absz1<- abs(x1d-x1s)
absz2<- abs(x2d-x2s)
absZ<- absz1+absz2

#store data
df<-data.frame(p1data,p2data,x1d,x2d,x1s,x2s,absz1,absz2,absZ)

Having run this code, lets talk about finding equilibrium and generating the visual we see featured above.

Finding Equlibrium

The primary way in which we will find a general equilibrium (an equilibrium in which both markets clear) in this model is through searching for minimum sum of absolute excess demands (in our code above this is absZ). Our finder is simply just looking for the row with this minimum value. That is:

#Find equilibrium
eqm<-subset(df,absZ==min(absZ))
View(eqm)

Which is:

Our Equilibrium Vector

Now this is not perfect in the sense that we will have the case where absZ=0, if we want to get a lower number we would have to choose a much more granular sequence of prices. This however comes at the cost of time in computation and the “curse of dimensionality”. Fortunately we can get sufficiently close where this issue is negligible.

Code for the Visual

The code for the visual is quite involved in terms of using ggplot(). Fortunately this is lots of copy and paste work with minor changes to the variables used. The one package which helps a lot here is the ggpubr package along with its ggarrage() function. This can all be seen below:

#Plot it in ggplot
p1x1<-ggplot(df,aes(y=df$p1data))+geom_point(aes(x=df$x1d),color="blue")+geom_point(aes(x=df$x1s),color="red")+xlab("Good 1")+ylab("Price 1")+theme_solarized()
p1x2<-ggplot(df,aes(y=df$p1data))+geom_point(aes(x=df$x2d),color="blue")+geom_point(aes(x=df$x2s),color="red")+xlab("Good 2")+ylab("Price 1")+theme_solarized()
p2x1<-ggplot(df,aes(y=df$p2data))+geom_point(aes(x=df$x1d),color="blue")+geom_point(aes(x=df$x1s),color="red")+xlab("Good 1")+ylab("Price 2")+theme_solarized()
p2x2<-ggplot(df,aes(y=df$p2data))+geom_point(aes(x=df$x2d),color="blue")+geom_point(aes(x=df$x2s),color="red")+xlab("Good 2")+ylab("Price 2")+theme_solarized()


ggarrange(p1x1,p1x2,p2x1,p2x2,ncol=2,nrow=2)

Conclusion

Coding general equilibrium models in R require baby steps. Using optimization functions and for loops have worked best for me so far. Some things which would be interesting to see is are economies with different preferences for consumers and different competitive environments faced by firms. I hope to work on more complex models and keep trying to improve my best practices. Let me know your thoughts!

Advertisement

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s

%d bloggers like this: