# Empirical project 11


## Getting started in Python


Head to the ‘Getting started in Python’ page for help and advice on setting up a Python session to work with. Remember, you can run any page from this book as a *notebook* by downloading the relevant file from this [repository](https://github.com/aeturrell/core_python) and running it on your own computer. Alternatively, you can run pages online in your browser over at [Binder](https://mybinder.org/v2/gh/aeturrell/core_python/HEAD).


## Preliminary settings


Let's import the packages we'll need and also configure the settings we want:


```python
import pandas as pd
import matplotlib as mpl
import matplotlib.pyplot as plt
import numpy as np
from pathlib import Path
import pingouin as pg
from lets_plot import *
from lets_plot.mapping import as_discrete


LetsPlot.setup_html(no_js=True)
```


### *Python walk-through 11.1* Importing data and recoding variables


Before importing data in Excel or .csv format, open it to ensure you understand the structure of the data and check if any additional options are required for the `pd.read_excel` function in order to import the data correctly. In this case the data is in a worksheet called ‘Data’, there are no missing values to worry about, and the first row contains the variable names. This format should be straightforward to import *but* this file is in the ‘Strict Open XML Spreadsheet (.xlsx)’ format rather than the ‘Excel Workbook (.xlsx)’ format, so loading it with `pd.read_excel` will fail. You'll need to re-save it in ‘Excel Workbook (.xlsx)’ format before importing the data using the `pd.read_excel` function.


```python
wtp = pd.read_excel(
    Path("data/doing-economics-datafile-working-in-excel-project-11.xlsx"),
    sheet_name="Data",
)
wtp.head(3)
```


<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }


```
.dataframe tbody tr th {
    vertical-align: top;
}


.dataframe thead th {
    text-align: right;
}
```


</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>id</th>
      <th>sex</th>
      <th>age</th>
      <th>education</th>
      <th>voc_training</th>
      <th>prob_a</th>
      <th>prob_b</th>
      <th>prob_c</th>
      <th>prob_d</th>
      <th>prob_e</th>
      <th>...</th>
      <th>cog_4</th>
      <th>cog_5</th>
      <th>cog_6</th>
      <th>liv_sit</th>
      <th>region</th>
      <th>occ_sit</th>
      <th>kids_nr</th>
      <th>member</th>
      <th>party</th>
      <th>hhnetinc</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>0</th>
      <td>2</td>
      <td>female</td>
      <td>50 - 59</td>
      <td>4</td>
      <td>3.0</td>
      <td>0</td>
      <td>0</td>
      <td>1</td>
      <td>0</td>
      <td>0</td>
      <td>...</td>
      <td>5</td>
      <td>4</td>
      <td>5</td>
      <td>not living with partner in HH</td>
      <td>Hessen</td>
      <td>Full-time employment</td>
      <td>no children</td>
      <td>no</td>
      <td>SPD</td>
      <td>1500 bis unter 2000 Euro</td>
    </tr>
    <tr>
      <th>1</th>
      <td>6</td>
      <td>female</td>
      <td>40 - 49</td>
      <td>6</td>
      <td>5.0</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>0</td>
      <td>1</td>
      <td>...</td>
      <td>4</td>
      <td>3</td>
      <td>4</td>
      <td>married and living with spouse</td>
      <td>Hamburg</td>
      <td>Occasionally or irregularly employed</td>
      <td>one child</td>
      <td>no</td>
      <td>Die Linke</td>
      <td>2600 bis unter 3200 Euro</td>
    </tr>
    <tr>
      <th>2</th>
      <td>8</td>
      <td>male</td>
      <td>50 - 59</td>
      <td>6</td>
      <td>6.0</td>
      <td>0</td>
      <td>0</td>
      <td>1</td>
      <td>0</td>
      <td>0</td>
      <td>...</td>
      <td>3</td>
      <td>3</td>
      <td>2</td>
      <td>married and living with spouse</td>
      <td>Nordrhein-Westfalen</td>
      <td>Full-time employment</td>
      <td>two children</td>
      <td>no</td>
      <td>Piraten</td>
      <td>5000 bis unter 6000 Euro</td>
    </tr>
  </tbody>
</table>
<p>3 rows × 79 columns</p>
</div>


#### Reverse-code variables


The first task is to recode variables related to the respondents’ views on certain aspects of government behaviour and attitudes about global warming (`cog_2`, `cog_5`, `scepticism_6`, and `scepticism_7`). This coding makes the interpretation of high and low values consistent across all questions, since the survey questions do not have this consistency.


To recode all of these values (across several variables) in one go, we use dictionary mapping: that is, using the `map` function to convert specific values to new values.


```python
map_dict = {1: 5, 2: 4, 3: 3, 4: 2, 5: 1}
vars_of_interest = ["cog_2", "cog_5", "scepticism_6", "scepticism_7"]
for col in vars_of_interest:
    wtp[col] = wtp[col].map(map_dict)
```


#### Create new variables containing WTP amounts


Although we could employ the same dictionary technique to recode the values for the minimum and maximum willingness to pay variables, we could use the `pd.merge` function instead. This function allows us to combine two dataframes via values given in a particular variable.


The `pd.merge` function is an essential tool for doing applied economics in Python. The `how` and `on` methods for this function are important to understand, so make sure to use the `help` function (run `help(pd.merge)`) to figure out what the different options for these methods are. You can find further information at the excellent **pandas** [documentation](https://pandas.pydata.org/docs/user_guide/merging.html#merging), via this [post](https://towardsdatascience.com/why-and-how-to-use-merge-with-pandas-in-python-548600f7e738), or in the online book [Coding for Economists](https://aeturrell.github.io/coding-for-economists/intro.html) (search for `pd.merge`).


We start by creating a new dataframe (`category_amount`) that has two variables: the original category value and the corresponding new euro amount. We then apply the `pd.merge` function to the `wtp` dataframe and the new dataframe, specifying the variables that link the data in each dataframe together. We're going to match on the `WTP_plmin` and, separately, the `WTP_plmax` options, to put in new columns for the min WTP and max WTP (`WTP_plmin_euro` and `WTP_plmax_euro` respectively).


```python
# vector containing the euro amounts
wtp_euro_levels = [48, 72, 84, 108, 156, 192, 252, 324, 432, 540, 720, 960, 1200, 1440]


# create a dataframe from this vector
category_amount = pd.DataFrame({"original": range(1, 15), "new": wtp_euro_levels})


# creating a new column for min WTP
wtp = pd.merge(
    category_amount, wtp, how="right", left_on="original", right_on="WTP_plmin"
).rename(columns={"new": "WTP_plmin_euro"})


# creating a new column for max WTP
wtp = pd.merge(
    category_amount, wtp, how="right", left_on="original", right_on="WTP_plmax"
).rename(columns={"new": "WTP_plmax_euro"})
```


Note that after the merging process you will find the column on which the two dataframes were merged (for example, `"original"` and `"WTP_plain"`) twice in the dataset. You may wish to delete one of them to avoid any subsequent confusion!


### *Python walk-through 11.2* Creating indices from multiple columns


We can create all of the required indices in three steps using column-aggregation operations; here the `.mean`. In each step we use the relevant function with `axis=1`, that is, aggregation over columns. We then insert the single, created column back into the dataframe. Our new columns are called `climate`, `gov_intervention`, and `pro_environment`. 


```python
wtp["climate"] = wtp[["scepticism_2", "scepticism_6", "scepticism_7"]].mean(axis=1)
wtp["gov_intervention"] = wtp[
    ["cog_1", "cog_2", "cog_3", "cog_4", "cog_5", "cog_6"]
].mean(axis=1)
wtp["pro_environment"] = wtp[["PN_1", "PN_2", "PN_3", "PN_4", "PN_6", "PN_7"]].mean(
    axis=1
)
```


### *Python walk-through 11.3* Calculating correlation coefficients


#### Calculate correlation coefficients and Cronbach's alpha


We showed how to calculate correlation coefficients in Python walk-through 10.1. In this case, since there are no missing values, we can use the `.corr` method without any additional options. We will round all values to 2 decimal places. 


For the questions on climate change:


```python
wtp[["scepticism_2", "scepticism_6", "scepticism_7"]].corr().round(2)
```


<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }


```
.dataframe tbody tr th {
    vertical-align: top;
}


.dataframe thead th {
    text-align: right;
}
```


</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>scepticism_2</th>
      <th>scepticism_6</th>
      <th>scepticism_7</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>scepticism_2</th>
      <td>1.00</td>
      <td>0.39</td>
      <td>0.42</td>
    </tr>
    <tr>
      <th>scepticism_6</th>
      <td>0.39</td>
      <td>1.00</td>
      <td>0.46</td>
    </tr>
    <tr>
      <th>scepticism_7</th>
      <td>0.42</td>
      <td>0.46</td>
      <td>1.00</td>
    </tr>
  </tbody>
</table>
</div>


For the questions on government behaviour:


```python
wtp[["cog_1", "cog_2", "cog_3", "cog_4", "cog_5", "cog_6"]].corr().round(2)
```


<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }


```
.dataframe tbody tr th {
    vertical-align: top;
}


.dataframe thead th {
    text-align: right;
}
```


</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>cog_1</th>
      <th>cog_2</th>
      <th>cog_3</th>
      <th>cog_4</th>
      <th>cog_5</th>
      <th>cog_6</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>cog_1</th>
      <td>1.00</td>
      <td>0.25</td>
      <td>0.32</td>
      <td>0.68</td>
      <td>0.29</td>
      <td>0.41</td>
    </tr>
    <tr>
      <th>cog_2</th>
      <td>0.25</td>
      <td>1.00</td>
      <td>0.12</td>
      <td>0.28</td>
      <td>0.41</td>
      <td>0.08</td>
    </tr>
    <tr>
      <th>cog_3</th>
      <td>0.32</td>
      <td>0.12</td>
      <td>1.00</td>
      <td>0.33</td>
      <td>0.02</td>
      <td>0.31</td>
    </tr>
    <tr>
      <th>cog_4</th>
      <td>0.68</td>
      <td>0.28</td>
      <td>0.33</td>
      <td>1.00</td>
      <td>0.27</td>
      <td>0.46</td>
    </tr>
    <tr>
      <th>cog_5</th>
      <td>0.29</td>
      <td>0.41</td>
      <td>0.02</td>
      <td>0.27</td>
      <td>1.00</td>
      <td>0.10</td>
    </tr>
    <tr>
      <th>cog_6</th>
      <td>0.41</td>
      <td>0.08</td>
      <td>0.31</td>
      <td>0.46</td>
      <td>0.10</td>
      <td>1.00</td>
    </tr>
  </tbody>
</table>
</div>


And, finally, the questions on personal behaviour:


```python
wtp[["PN_1", "PN_2", "PN_3", "PN_4", "PN_6", "PN_7"]].corr().round(2)
```


<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }


```
.dataframe tbody tr th {
    vertical-align: top;
}


.dataframe thead th {
    text-align: right;
}
```


</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>PN_1</th>
      <th>PN_2</th>
      <th>PN_3</th>
      <th>PN_4</th>
      <th>PN_6</th>
      <th>PN_7</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>PN_1</th>
      <td>1.00</td>
      <td>0.48</td>
      <td>0.43</td>
      <td>0.42</td>
      <td>0.41</td>
      <td>0.46</td>
    </tr>
    <tr>
      <th>PN_2</th>
      <td>0.48</td>
      <td>1.00</td>
      <td>0.63</td>
      <td>0.44</td>
      <td>0.50</td>
      <td>0.65</td>
    </tr>
    <tr>
      <th>PN_3</th>
      <td>0.43</td>
      <td>0.63</td>
      <td>1.00</td>
      <td>0.46</td>
      <td>0.52</td>
      <td>0.59</td>
    </tr>
    <tr>
      <th>PN_4</th>
      <td>0.42</td>
      <td>0.44</td>
      <td>0.46</td>
      <td>1.00</td>
      <td>0.57</td>
      <td>0.39</td>
    </tr>
    <tr>
      <th>PN_6</th>
      <td>0.41</td>
      <td>0.50</td>
      <td>0.52</td>
      <td>0.57</td>
      <td>1.00</td>
      <td>0.46</td>
    </tr>
    <tr>
      <th>PN_7</th>
      <td>0.46</td>
      <td>0.65</td>
      <td>0.59</td>
      <td>0.39</td>
      <td>0.46</td>
      <td>1.00</td>
    </tr>
  </tbody>
</table>
</div>


#### Calculate Cronbach's alpha


We can compute Cronbach’s alpha using the `cronbach_alpha` function from the **pingouin** package. Let's calculate it for these three sets of data:


```python
pg.cronbach_alpha(data=wtp[["scepticism_2", "scepticism_6", "scepticism_7"]])
```


```
(0.684752824338845, array([0.656, 0.711]))
```


```python
pg.cronbach_alpha(data=wtp[["cog_1", "cog_2", "cog_3", "cog_4", "cog_5", "cog_6"]])
```


```
(0.7133064708899756, array([0.69 , 0.735]))
```


```python
pg.cronbach_alpha(data=wtp[["PN_1", "PN_2", "PN_3", "PN_4", "PN_6", "PN_7"]])
```


```
(0.8521220890874502, array([0.84 , 0.863]))
```


### *Python walk-through 11.4* Using loops to obtain summary statistics


The two different formats (DC and TWPL) are recorded in the variable `abst_format`, and take the values `ref` and `ladder` respectively. We will store all the variables of interest into a list called `variables`, then use a ‘for’ loop to calculate summary statistics for each variable and present it in a table (using `pd.crosstab`).


```python
variables = ["sex", "age", "kids_nr", "hhnetinc", "member", "education"]


for var in variables:
    print(pd.crosstab(wtp[var], wtp["abst_format"], normalize="columns").round(2))
    print("————\n")
```


```
abst_format  ladder   ref
sex                      
female         0.52  0.52
male           0.48  0.48
————


abst_format  ladder   ref
age                      
18 - 24        0.09  0.10
25 - 29        0.08  0.09
30 - 39        0.18  0.17
40 - 49        0.22  0.23
50 - 59        0.24  0.24
60 - 69        0.18  0.18
————


abst_format            ladder   ref
kids_nr                            
four or more children    0.01  0.01
no children              0.65  0.66
one child                0.20  0.18
three children           0.03  0.03
two children             0.11  0.12
————


abst_format               ladder   ref
hhnetinc                              
1100 bis unter 1500 Euro    0.14  0.13
1500 bis unter 2000 Euro    0.15  0.15
2000 bis unter 2600 Euro    0.11  0.15
2600 bis unter 3200 Euro    0.11  0.11
3200 bis unter 4000 Euro    0.11  0.08
4000 bis unter 5000 Euro    0.05  0.05
500 bis unter 1100 Euro     0.13  0.14
5000 bis unter 6000 Euro    0.03  0.02
6000 bis unter 7500 Euro    0.01  0.00
7500 und mehr               0.00  0.00
bis unter 500 Euro          0.03  0.04
do not want to answer       0.12  0.13
————


abst_format  ladder   ref
member                   
no             0.92  0.91
yes            0.08  0.09
————


abst_format  ladder   ref
education                
1              0.01  0.01
2              0.02  0.02
3              0.34  0.33
4              0.26  0.27
5              0.07  0.07
6              0.29  0.30
————
```


The output above gives the required tables, but it is not easy to read. You may want to tidy up the results, for example by translating the household income amounts (`hhnetinc`) from German to English and reordering the options from smallest to largest.


### *Python walk-through 11.5* Calculating summary statistics


The `agg` function can provide multiple statistics for a number of variables in one command. You will need to provide a list of the variables you want to summarise (after your `groupby`), then use the `agg` option to specify the summary statistics you need. Here, we need the mean, standard deviation, min, and max for the variables `climate`, `gov_intervention`, and `pro_environment`. Finally, `.stack` puts the data into a longer (and here more readable) format.


```python
(
    wtp.groupby(["abst_format"])[["climate", "gov_intervention", "pro_environment"]]
    .agg(["mean", "std", "min", "max"])
    .stack()
)
```


<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }


```
.dataframe tbody tr th {
    vertical-align: top;
}


.dataframe thead th {
    text-align: right;
}
```


</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th></th>
      <th>climate</th>
      <th>gov_intervention</th>
      <th>pro_environment</th>
    </tr>
    <tr>
      <th>abst_format</th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th rowspan="4" valign="top">ladder</th>
      <th>mean</th>
      <td>2.285903</td>
      <td>3.145586</td>
      <td>3.030962</td>
    </tr>
    <tr>
      <th>std</th>
      <td>0.836167</td>
      <td>0.698730</td>
      <td>0.792378</td>
    </tr>
    <tr>
      <th>min</th>
      <td>1.000000</td>
      <td>1.000000</td>
      <td>1.000000</td>
    </tr>
    <tr>
      <th>max</th>
      <td>5.000000</td>
      <td>5.000000</td>
      <td>5.000000</td>
    </tr>
    <tr>
      <th rowspan="4" valign="top">ref</th>
      <th>mean</th>
      <td>2.373426</td>
      <td>3.187210</td>
      <td>3.012591</td>
    </tr>
    <tr>
      <th>std</th>
      <td>0.847227</td>
      <td>0.658170</td>
      <td>0.819205</td>
    </tr>
    <tr>
      <th>min</th>
      <td>1.000000</td>
      <td>1.000000</td>
      <td>1.000000</td>
    </tr>
    <tr>
      <th>max</th>
      <td>5.000000</td>
      <td>5.000000</td>
      <td>5.000000</td>
    </tr>
  </tbody>
</table>
</div>


### *Python walk-through 11.6* Summarising willingness to pay variables


#### Create column charts for minimum and maximum WTP


Before we can plot a column chart, we need to compute frequencies (number of observations) for each value of the willingness to pay (1–14). We do this separately for the minimum and maximum willingness to pay.


In each case we select the relevant variable and remove any observations with missing values using the `.dropna()` method. We can then separate the data by level (WTP amount) of the `WTP_plmin_euro` or `WTP_plmax_euro` variables (using `groupby`), then obtain a frequency count using the `todo` function.


Once we have the frequency count stored as a dataframe, we can plot the column charts.


For the minimum willingness to pay:


```python
df_plmin = (
    wtp[["WTP_plmin_euro"]]
    .dropna()
    .astype("int")
    .astype("category")
    .value_counts()
    .sort_values()
    .reset_index()
)
```


```python
(
    ggplot(
        df_plmin.sort_values(by="WTP_plmin_euro"),
        aes(x=as_discrete("WTP_plmin_euro"), y="count"),
    )
    + geom_bar(stat="identity")
    + labs(
        x="Minimum WTP (euros)",
        y="Frequency",
    )
)
```


![png](empirical_project_11_files/empirical_project_11_WT6_1.png)


**Figure 11.4** Minimum WTP (euros).


Let's now do the same for the maximum willingness to pay:


```python
df_plmax = (
    wtp[["WTP_plmax_euro"]]
    .dropna()
    .astype("int")
    .astype("category")
    .value_counts()
    .sort_values()
    .reset_index()
)
```


```python
(
    ggplot(df_plmax, aes(x=as_discrete("WTP_plmax_euro"), y="count"))
    + geom_bar(stat="identity")
    + labs(x="Maximum WTP (euros)", y="Frequency")
)
```


![png](empirical_project_11_files/empirical_project_11_WT6_2.png)


**Figure 11.5** Maximum WTP (euros).


#### Calculate average WTP for each individual


We can use the `mean` function to obtain the average of the minimum and maximum willingness to pay (combining the two columns at each row using the `axis=1` keyword argument). We store the averages as the variable `wtp_average`. 


```python
wtp["wtp_average"] = wtp[["WTP_plmin_euro", "WTP_plmax_euro"]].mean(axis=1)
```


#### Calculate mean and median WTP across individuals


The mean and median of this average value can be obtained using the `mean` and `median` functions. Note that invalid entries, such as `NaN`, are omitted by default. Also, since there's only a single remaining dimension over which to take the mean and median, there’s no need to specify `axis=0`, which is the default in any case.


```python
wtp["wtp_average"].mean()
```


```
268.5344827586207
```


```python
wtp["wtp_average"].median()
```


```
132.0
```


#### Calculate correlation coefficients


We showed how to obtain a matrix of correlation coefficients for a number of variables in Python walk-through 8.8. We use the same process here, storing the coefficients in an object called `M_corr`.


```python
M_corr = (
    wtp.assign(gender=lambda x: np.where(x["sex"] == "female", 0, 1))
    .loc[
        :,
        [
            "wtp_average",
            "education",
            "gender",
            "climate",
            "gov_intervention",
            "pro_environment",
        ],
    ]
    .corr()
)


M_corr["wtp_average"].round(3)
```


```
wtp_average         1.000
education           0.138
gender              0.037
climate            -0.145
gov_intervention   -0.188
pro_environment     0.188
Name: wtp_average, dtype: float64
```


### *Python walk-through 11.7* Summarising dichotomous choice (DC) variables


#### Create a frequency table for `DC_ref_outcome`


We can group by `costs` and `DC_ref_outcome` to obtain the number of observations for each combination of amount and vote response. We can also recode the voting options as shorter phrases: ‘Yes’, ‘No’, and ‘Abstain’.


```python
recoding_dict = {
    "do not support referendum and no pay": "No",
    "support referendum and pay": "Yes",
    "would not vote": "Abstain",
}


wtp_dc = (
    wtp.dropna(subset=["costs", "DC_ref_outcome"])
    .assign(DC_ref_outcome=lambda x: x["DC_ref_outcome"].map(recoding_dict))
    .groupby(["costs", "DC_ref_outcome"])["id"]
    .count()
    .unstack()
)


wtp_dc
```


<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }


```
.dataframe tbody tr th {
    vertical-align: top;
}


.dataframe thead th {
    text-align: right;
}
```


</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th>DC_ref_outcome</th>
      <th>Abstain</th>
      <th>No</th>
      <th>Yes</th>
    </tr>
    <tr>
      <th>costs</th>
      <th></th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>48</th>
      <td>12</td>
      <td>21</td>
      <td>32</td>
    </tr>
    <tr>
      <th>72</th>
      <td>11</td>
      <td>30</td>
      <td>40</td>
    </tr>
    <tr>
      <th>84</th>
      <td>12</td>
      <td>24</td>
      <td>45</td>
    </tr>
    <tr>
      <th>108</th>
      <td>7</td>
      <td>35</td>
      <td>31</td>
    </tr>
    <tr>
      <th>156</th>
      <td>13</td>
      <td>31</td>
      <td>40</td>
    </tr>
    <tr>
      <th>192</th>
      <td>11</td>
      <td>25</td>
      <td>25</td>
    </tr>
    <tr>
      <th>252</th>
      <td>9</td>
      <td>32</td>
      <td>28</td>
    </tr>
    <tr>
      <th>324</th>
      <td>16</td>
      <td>41</td>
      <td>27</td>
    </tr>
    <tr>
      <th>432</th>
      <td>11</td>
      <td>35</td>
      <td>29</td>
    </tr>
    <tr>
      <th>540</th>
      <td>9</td>
      <td>31</td>
      <td>22</td>
    </tr>
    <tr>
      <th>720</th>
      <td>12</td>
      <td>39</td>
      <td>13</td>
    </tr>
    <tr>
      <th>960</th>
      <td>14</td>
      <td>28</td>
      <td>15</td>
    </tr>
    <tr>
      <th>1200</th>
      <td>11</td>
      <td>42</td>
      <td>21</td>
    </tr>
    <tr>
      <th>1440</th>
      <td>19</td>
      <td>42</td>
      <td>15</td>
    </tr>
  </tbody>
</table>
</div>


#### Add column showing the proportion voting yes or no


We can extend the table from Question 2(a) to include the proportion voting yes or no (to obtain percentages, multiply the values by 100). We *chain* methods in the below; using `assign` to create new columns and `round` to round the values in the float columns.


```python
wtp_dc = wtp_dc.assign(
    total=lambda x: x["Abstain"] + x["No"] + x["Yes"],
    prop_no=lambda x: (x["Abstain"] + x["No"]) / x["total"],
    prop_yes=lambda x: x["Yes"] / x["total"],
).round(2)


wtp_dc
```


<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }


```
.dataframe tbody tr th {
    vertical-align: top;
}


.dataframe thead th {
    text-align: right;
}
```


</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th>DC_ref_outcome</th>
      <th>Abstain</th>
      <th>No</th>
      <th>Yes</th>
      <th>total</th>
      <th>prop_no</th>
      <th>prop_yes</th>
    </tr>
    <tr>
      <th>costs</th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>48</th>
      <td>12</td>
      <td>21</td>
      <td>32</td>
      <td>65</td>
      <td>0.51</td>
      <td>0.49</td>
    </tr>
    <tr>
      <th>72</th>
      <td>11</td>
      <td>30</td>
      <td>40</td>
      <td>81</td>
      <td>0.51</td>
      <td>0.49</td>
    </tr>
    <tr>
      <th>84</th>
      <td>12</td>
      <td>24</td>
      <td>45</td>
      <td>81</td>
      <td>0.44</td>
      <td>0.56</td>
    </tr>
    <tr>
      <th>108</th>
      <td>7</td>
      <td>35</td>
      <td>31</td>
      <td>73</td>
      <td>0.58</td>
      <td>0.42</td>
    </tr>
    <tr>
      <th>156</th>
      <td>13</td>
      <td>31</td>
      <td>40</td>
      <td>84</td>
      <td>0.52</td>
      <td>0.48</td>
    </tr>
    <tr>
      <th>192</th>
      <td>11</td>
      <td>25</td>
      <td>25</td>
      <td>61</td>
      <td>0.59</td>
      <td>0.41</td>
    </tr>
    <tr>
      <th>252</th>
      <td>9</td>
      <td>32</td>
      <td>28</td>
      <td>69</td>
      <td>0.59</td>
      <td>0.41</td>
    </tr>
    <tr>
      <th>324</th>
      <td>16</td>
      <td>41</td>
      <td>27</td>
      <td>84</td>
      <td>0.68</td>
      <td>0.32</td>
    </tr>
    <tr>
      <th>432</th>
      <td>11</td>
      <td>35</td>
      <td>29</td>
      <td>75</td>
      <td>0.61</td>
      <td>0.39</td>
    </tr>
    <tr>
      <th>540</th>
      <td>9</td>
      <td>31</td>
      <td>22</td>
      <td>62</td>
      <td>0.65</td>
      <td>0.35</td>
    </tr>
    <tr>
      <th>720</th>
      <td>12</td>
      <td>39</td>
      <td>13</td>
      <td>64</td>
      <td>0.80</td>
      <td>0.20</td>
    </tr>
    <tr>
      <th>960</th>
      <td>14</td>
      <td>28</td>
      <td>15</td>
      <td>57</td>
      <td>0.74</td>
      <td>0.26</td>
    </tr>
    <tr>
      <th>1200</th>
      <td>11</td>
      <td>42</td>
      <td>21</td>
      <td>74</td>
      <td>0.72</td>
      <td>0.28</td>
    </tr>
    <tr>
      <th>1440</th>
      <td>19</td>
      <td>42</td>
      <td>15</td>
      <td>76</td>
      <td>0.80</td>
      <td>0.20</td>
    </tr>
  </tbody>
</table>
</div>


#### Make a line chart of WTP


Using the dataframe generated for Questions 2(a) and (b) (`wtp_dc`), we can use **lets-plot** to plot the **demand curve** as a scatterplot with connected points.


```python
(
    ggplot(wtp_dc.reset_index(), aes(x="costs", y="prop_yes"))
    + geom_line(size=1)
    + geom_point(size=4)
    + labs(
        x="Amount (euros)",
        y="Proportion voting 'Yes'",
    )
)
```


![png](empirical_project_11_files/empirical_project_11_WT7_1.png)


**Figure 11.6** Demand curve (in euros), DC method.


#### Calculate new proportions and add them to the table and chart


It is straightforward to calculate new proportions and add them to the existing dataframe:


```python
wtp_dc = wtp_dc.assign(
    total_ex=lambda x: x["No"] + x["Yes"],
    prop_no_ex=lambda x: (x["No"]) / x["total_ex"],
    prop_yes_ex=lambda x: x["Yes"] / x["total_ex"],
).round(2)
```


Now we're going to plot them too. Because **lets-plot** expects long-format (aka ‘tidy’) data, we need to re-orient the dataframe. We'll put the results of this long format data in a new dataframe called `demand_curve`. We'll clean this dataframe up a bit too, by renaming the ref outcome variable to `vote` and the entries in that variable to more relevant names for the chart. Finally, we'll plot the chart using lines and dots as before.


```python
demand_curve = (
    pd.melt(
        wtp_dc.reset_index(), id_vars=["costs"], value_vars=["prop_yes", "prop_yes_ex"]
    )
    .rename(columns={"DC_ref_outcome": "Vote"})
    .assign(
        Vote=lambda x: x["Vote"].map(
            {"prop_yes": "Counted as no", "prop_yes_ex": "Excluded"}
        )
    )
)
```


```python
(
    ggplot(demand_curve, aes(x="costs", y="value", color="vote"))
    + geom_line(size=1)
    + geom_point(size=3)
    + labs(
        x="Amount (euros)",
        y="Proportion voting 'Yes'",
    )
    + ggsize(900, 400)
)
```


** Figure here ** 
![png](empirical_project_11_files/empirical_project_11_WT7_2.png)


**Figure 11.7** Demand curve from DC respondents, under different treatments for 'Abstain' responses.


### *Python walk-through 11.8* Calculating confidence intervals for differences in means


#### Calculate the difference in means, standard deviations, and number of observations


We first create two vectors that will contain the `wtp` values for each of the two question methods. For the DC format, willingness to pay is recorded in the `costs` variable, so we select all observations where the `DC_ref_outcome` variable indicates the individual voted `yes` and drop any missing observations. For the TWPL format we use the `wtp_average` variable that we created in Python walk-through 11.6.


```python
dc_wtp = wtp.loc[
    wtp["DC_ref_outcome"] == "support referendum and pay", "costs"
].dropna()


dc_wtp.agg(["mean", "std", "count"]).round(2)
```


```
mean     348.19
std      378.65
count    383.00
Name: costs, dtype: float64
```


```python
twpl_wt = wtp.loc[~wtp["wtp_average"].isna(), "wtp_average"].dropna()


twpl_wt.agg(["mean", "std", "count"]).round(2)
```


```
mean     268.53
std      287.71
count    348.00
Name: wtp_average, dtype: float64
```


#### Calculate 95% confidence intervals


Using the `ttest` function from the **pingouin** package to obtain 95% confidence intervals was covered in Python walk-throughs 8.10 and 10.6. As we have already separated the data for the two different question formats in Question 3(a), we can obtain the confidence interval directly.


```python
pg.ttest(dc_wtp, twpl_wt, confidence=0.05)
```


<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }


```
.dataframe tbody tr th {
    vertical-align: top;
}


.dataframe thead th {
    text-align: right;
}
```


</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>T</th>
      <th>dof</th>
      <th>alternative</th>
      <th>p-val</th>
      <th>CI5%</th>
      <th>cohen-d</th>
      <th>BF10</th>
      <th>power</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>T-test</th>
      <td>3.219251</td>
      <td>707.306861</td>
      <td>two-sided</td>
      <td>0.001344</td>
      <td>[78.10141039468779, 81.20560320034237]</td>
      <td>0.235367</td>
      <td>12.893</td>
      <td>0.887627</td>
    </tr>
  </tbody>
</table>
</div>


#### Calculate median WTP for the DC format


In Python walk-through 11.6 we obtained the median WTP for the TWPL format. We now obtain the WTP using the DC format:


```python
dc_wtp.median()
```


```
192.0
```