# The WTP for marine protected areas: Replicating Karlõševa et al. (2016)

In the spirit of replicating studies, in this blog post I am replicating the findings from Karlõševa et al. (2016) and estimating the willingness to pay (WTP) to establish a marine protected area.

The aim of Karlõševa et al. (2016) is to investigate preferences for offshore wind farms and marine protected areas in the coast of Estonia, which are two policy options currently being considered by authorities. The authors postulate that public opposes offshore wind farms when compared with the establishment of marine protected areas.

Karlõševa et al. (2016) conduct a choice experiment, wherein individuals are asked to choose one out of three alternatives: whether to establish a wind farm, an “eco-friendly” wind farm, or a marine protected area in two study areas in the coast of Estonia. They conducted a survey to 800 individuals thorugh an Internet panel back in 2013. The authors then estimate a latent class mixed logit model, and estimate the WTP (in EUR) per year and per household. To establish a marine protected area the mean WTP is 29.13 EUR, -10.47 EUR for a wind farm, and 25.46 EUR for an eco-friendly wind farm.

Instead of estimating the latent class mixed logit model like Karlõševa et al. (2016) did, I will estimate a simpler conditional logit model using their data.

The data is available here under “Research Data”. To replicate these findings, you have to download the file called  “data_EE1_01.mat”.

This file is a matlab file. In order to import it to R, I use the R.matlab package and the readMat function in this package.

install.packages("R.matlab")
library(R.matlab)
DATA <- readMat("DIRECTORY/data_EE1_01.mat")

Unfortunately the readMat function imports the file as a list, rather than a dataframe. So it looks like the picture in the left, rather than the one on the right, which is what I want.

To get the dataframe as shown by the picture in the left, I will run the following command to convert the list DATA into a dataframe object, which I call Data.df:

n <- length(DATA[[1]])
Data.df <- structure(DATA, row.names = c(NA, -n), class = "data.frame")

Once I could get a look at the data, I identified some of the variables we need to estimate the conditional logit model:

• “ans1” is the choice made by the respondent, that is, 1 is if the alternative is chosen, and zero if it is not chosen
• “alt” refers to the alternative 1, 2 or 3 presented in each choice card, which is denoted by “cs3”
• “id” is the individual identifier
• “as” refers to the scenario in the Apollo shoal study area
• “ws” refers to the scenario in the Western shoal
• “c” is the cost vector, ranging from 0 to 20

I use the gmnl package to estimate the conditional logit model. In this previous blog post, I also used the gmnl package to handle this triary data. “as” and “ws” are our variables of interest. I am assuming that when these take the value 1, the respondents faced the marine protected area scenario; when it takes the value 2, the respondents faced the offshore wind scenario; and when it takes the value 3, the respondents faced the eco-friendly offshore wind scenario.

The higher the cost “c”, the less likely it will be that the alternative will be chosen. “alt1” captures the average utility associated with the status quo alternative, that is not establishing a marine protected area nor offshore wind farm. I had to create the alt1 variable (see below).

As usual, I first call the gmnl package. I then create a mlogit object called df.data, by declaring what the id.var, choice, and alt.var are.

library(gmnl)
df.data<-mlogit::mlogit.data(Data.df,
id.var="id",
choice="ans1",
shape="long",
alt.var = "alt",
sep="")

df.data$alt1 <- 0 df.data$alt1[df.data\$alt==1] <- 1
model <- gmnl(ans1 ~ c + factor(ws) + factor(as) | 0, data=df.data)
summary(model)

The function to estimate the conditional logit model is gmnl, and its syntax ressembles the lm function. This is the resulting output:

Coefficients:
Estimate Std. Error  z-value  Pr(>|z|)
c             -0.0320756  0.0024913 -12.8749 < 2.2e-16 ***
1:(intercept)  0.7886318  0.0575478  13.7039 < 2.2e-16 ***
factor(ws)1    0.3360950  0.0523005   6.4262 1.308e-10 ***
factor(ws)2   -0.0138831  0.0518240  -0.2679  0.788785
factor(ws)3    0.1622532  0.0501342   3.2364  0.001211 **
factor(as)1    0.1932240  0.0474185   4.0749 4.604e-05 ***
factor(as)2   -0.4957064  0.0527695  -9.3938 < 2.2e-16 ***
factor(as)3    0.2966563  0.0480431   6.1748 6.625e-10 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Optimization of log-likelihood by Newton-Raphson maximisation
Log Likelihood: -9301.6
Number of observations: 9600

Estimating the WTP to establish a marine protected area is relatively straighforward:

$WTP_M = \frac{\beta_M}{\beta_C}$

The coefficient associated with establishing a marine protected area is $beta_M$=0.336 for the Western shoal and $beta_M$=0.193 in the Apollo shoal. $beta_C$=-0.032 is the coefficient associated with the cost variable. I can interpret these as meaning that utility increases in the scenarios where the marine protected area is established, and decreases with increases in the cost of the policy.

Accordingly, the WTP to establish a marine protected area for the Western and Apollo shoals are:

wtp_ws<--coef(model)[3]/coef(model)[1]
wtp_as<--coef(model)[6]/coef(model)[1]
> wtp_ws
factor(ws)1
10.47823
> wtp_as
factor(as)1
6.024026

Hence, I conclude that according to the conditional logit model, the WTPs to establish a marine protected area in the Western and Apollo shoals is 10.5 and 6.02 EUR per year and per household, respectively.

Is the WTP to establish a MPA in Western shoals statistically higher than in the Apollo shoal? One could find out by comparing the confidence intervals of these WTP estimates. To get the confidence intervals, I need to calculate the standard deviation of each WTP. To do so, I have several alternatives but for simplicity sake I will use the delta method. The function is in the gmnl package, namely the msm::deltamethod function:

se_wtp_ws<-msm::deltamethod(~x3/x1, coef(model),vcov(model))
se_wtp_as<-msm::deltamethod(~x6/x1, coef(model),vcov(model))

If I have the standard deviation, I can calculate the 95% confidence interval using the critical t-value of 1.96:

> c(wtp_ws-1.96*se_wtp_ws, wtp_ws , wtp_ws+1.96*se_wtp_ws)
factor(ws)1 factor(ws)1 factor(ws)1
7.012324   10.478228   13.944132
> c(wtp_as-1.96*se_wtp_as, wtp_as , wtp_as+1.96*se_wtp_as)
factor(as)1 factor(as)1 factor(as)1
2.976386    6.024026    9.071666

The WTP to establish a marine protected area in Western shoal falls within the 95% confidence interval of [7.01 , 13.94]. This WTP for Apollo shoal falls within the 95% confidence interval of [2.97 , 9.07]. Since the confidence intervals overlap, I cannot conclude that the two WTP are statistically different.

Are there other alternatives to compare the different WTP estimates? That question will be the focus of future blogs.

References:

Karlõševa, A., Nõmmann, S., Nõmmann, T., Urbel-Piirsalu, E., Budziński, W., Czajkowski, M., & Hanley, N. (2016). Marine trade-offs: Comparing the benefits of off-shore wind farms and marine protected areas. Energy Economics55, 127-134.

Karlõševa, Aljona (2016), “Data for: Marine trade-offs: Comparing the benefits of off-shore wind farms and marine protected areas ”, Mendeley Data, v1http://dx.doi.org/10.17632/j4zwwz8gh7.1