CViM / src / plots / script.jl
script.jl
Raw
using Pkg
Pkg.activate("src/plots")
Pkg.instantiate()

# Load packages & code libraries
using DataFrames
using CSV
using Pipe
using Statistics
using CairoMakie
using QuantEcon
using EasyFit
using LaTeXStrings

include("utils.jl")

function load_data()
    adf = DataFrame()
    mdf = DataFrame()

    for scenario in SCENARIOS, shock in SHOCKS
        append!(adf, CSV.File("data/shock=$(shock)/$(scenario)/adf.csv"); promote = true)
        append!(mdf, CSV.File("data/shock=$(shock)/$(scenario)/mdf.csv"); promote = true)
    end
    
    return adf, mdf
end

df, mdf = load_data();

const vars_ib = [:lending_facility, :deposit_facility, :liq_pref, :ON_interbank, :Term_interbank,
    :ib_on_demand, :ib_term_demand, :flow, :loans]
const vars_hh = [:networth, :consumption]
const vars_firms = [:networth, :consumption, :GDP]
const number_of_households = mean(mdf.n_hh)
const number_of_firms =  mean(mdf.n_f)

cd(mkpath("img/pdf")) do
    # Figure 3
    ## Figure 3a
    function plot_figure_3a_C6(df, vars, ylabels)
        fig = Figure(size = (800, 400))
        gdf = groupby(df, :scenario)

        for i in 1:length(vars)
            ax = fig[1, i] = Axis(fig, title = ylabels[i])
            for j in 1:length(gdf)
                # Apply HP filter to the selected variable
                cycle, trend = hp_filter(gdf[j][!, vars[i]], 1600)

                # Plot the main trend line
                lines!(trend; linewidth = 2, label = only(unique(gdf[j].scenario)), color =  j == 1 ? :black : Makie.wong_colors()[j])

                # Add standard deviation lines
                standard_deviation_bands!(cycle, trend, j == 1 ? :black : Makie.wong_colors()[j])
                
                # Set ticks x-axys
                ax.xticks = 0:200:800
                fig.content[1].ylabel = "Mean"
                fig.content[i].xlabel = "Steps"
                linkyaxes!(fig.content...)
                invisible_yaxis!(fig, i)
            end
            # Legend
        end
       fig[end+1,1:2] = Legend(
            fig,
            fig.content[1];
            tellheight = true,
            tellwidth = false,
            orientation = :horizontal
        )
        return fig
    end

    df_banks = generate_df(df, mdf; shock = "Missing", banks = true);
    p = plot_figure_3a_C6(filter(:status => x -> x == "borrower", df_banks), [:ON_interbank, :Term_interbank], ["ON volumes", "Term volumes"])
    save("Fig3a.pdf", p)

    ## Figure 3b
    function plot_figure_3b(df, vars_nominator, vars_denominator, ylabels)
        fig = Figure(size = (800, 400))
        gdf = groupby(df, :scenario)

        for i in 1:length(vars_nominator)
            ax = fig[1, i] = Axis(fig, title = ylabels[i])
            for j in 1:length(gdf)
                # Apply HP filter to the selected variable
                cycle, trend = hp_filter((1 .- (gdf[j][!, vars_nominator[i]] ./ gdf[j][!, vars_denominator[i]])) .* 100, 1600)

                # Plot the main trend line
                lines!(trend; linewidth = 2, label = only(unique(gdf[j].scenario)), color =  j == 1 ? :black : Makie.wong_colors()[j])

                # Add standard deviation lines
                standard_deviation_bands!(cycle, trend, j == 1 ? :black : Makie.wong_colors()[j])
                
                # Set ticks x-axys
                ax.xticks = 0:200:800
                fig.content[1].ylabel = "Rate (%)"
                fig.content[i].xlabel = "Steps"
                linkyaxes!(fig.content...)
                invisible_yaxis!(fig, i)
            end
            # Legend
        end
        fig[end+1,1:length(vars_nominator)] = Legend(
            fig,
            fig.content[1];
            tellheight = true,
            tellwidth = false,
            orientation = :horizontal
        )
        return fig
    end

    p = plot_figure_3b(filter(:status => ==("borrower"), df_banks), [:ON_interbank, :Term_interbank], [:ib_on_demand, :ib_term_demand], ["Overnight rationing", "Term rationing"])
    save("Fig3b.pdf", p)

    ## Figures 3c - 3d
    function plot_figures_3c_3d_C4_C5(dataframes, vars, ylabels)
        fig = Figure(size = (800, 400))

        for i in 1:length(dataframes)
            ax = fig[1, i] = Axis(fig, title = ylabels[i])
            # Take gdf for each dataframe
            gdf = @pipe dataframes[i] |>
                groupby(_, :scenario)
            for j in 1:length(gdf)
                # Apply HP filter to the selected variable
                cycle, trend = hp_filter(gdf[j][!, vars[1]], 1600)

                # Plot the main trend line
                lines!(trend; linewidth = 2, label = only(unique(gdf[j].scenario)), color =  j == 1 ? :black : Makie.wong_colors()[j])

                # Add standard deviation lines
                standard_deviation_bands!(cycle, trend, j == 1 ? :black : Makie.wong_colors()[j])

                # Set ticks x-axys
                ax.xticks = 0:200:800
                fig.content[1].ylabel = "Mean" 
                fig.content[i].xlabel = "Steps"
                linkyaxes!(fig.content...)
                invisible_yaxis!(fig, i)
            end
        end
        fig[end+1,1:length(dataframes)] = Legend(
            fig,
            fig.content[1];
            tellheight = true,
            tellwidth = false,
            orientation = :horizontal
        )
        return fig
    end

    p = plot_figures_3c_3d_C4_C5([filter(:status => ==("borrower"), df_banks), 
        filter(:status => ==("lender"), df_banks)], [:flow], ["Outflows", "Inflows"])
    save("Fig3c.pdf", p)

    p = plot_figures_3c_3d_C4_C5([filter(:status => ==("borrower"), df_banks), 
    filter(:status => ==("lender"), df_banks)], [:liq_pref], ["Borrowers liquidity preferences", "Lenders liquidity preferences"])
    save("Fig3d.pdf", p)

    # Figure 4
    function plots_by_group(df, vars)
        fig = Figure(size = (1200, 600))
        gdf = groupby(df, :scenario)

        for j in 1:length(gdf)
            axes = ((1:2,1), (1,2), (1,3), (2,2), (2,3))
            ax = fig[axes[j]...] = Axis(fig, title = SCENARIOS[j])

            for i in eachindex(VALUE_GROUPS)
                # Filter the data for the current status
                sdf = filter(r -> r.value == VALUE_GROUPS[i], gdf[j])

                # Apply HP filter to the selected variable
                cycle, trend = hp_filter(sdf[!, vars[1]], 1600)


                # Plot the main trend line
                lines!(trend; linewidth = 2, label = VALUE_GROUPS[i], color = Makie.wong_colors()[2:end][i])

                # Add standard deviation lines
                standard_deviation_bands!(cycle, trend, Makie.wong_colors()[2:end][i])

                # Set ticks x-axys
                ax.xticks = 0:200:800
                fig.content[1].ylabel = "Mean"
                fig.content[j].xlabel = j in [1, 4, 5] ? "Steps" : ""
                fig.content[j].xticksvisible = j in [2, 3] ? false : true  
                fig.content[j].xticklabelsvisible = j in [2, 3] ? false : true  
            end
        end
    configure_axes_links(fig, gdf)
    fig[end+1,1:3] = Legend(
            fig,
            fig.content[1];
            tellheight = true,
            tellwidth = false,
            orientation = :horizontal
        )
        return fig
    end

    df_banks_1 = generate_df(df, mdf; shock = "Missing", banks = true, value = true);
    p = plots_by_group(df_banks_1, [:ON_interbank])
    save("Fig4.pdf", p)

    # Figure 5
    function plot_figure_5(df, vars, labels)
        fig = Figure(size = (1200, 600))
        gdf = groupby(df, :scenario)

        # Iterate over each gdf
        for i in 1:length(gdf)
            # Define axes positions
            axes = ((1:2,1), (1,2), (1,3), (2,2), (2,3))
            ax = fig[axes[i]...] = Axis(fig, title = SCENARIOS[i])
            
            for j in 1:length(vars)
                # Apply HP filter to the selected variable
                cycle, trend = hp_filter(gdf[i][!, vars[j]], 1600)

                # Plot the main trend line with labels
                lines!(trend; label = labels[j], linewidth = 2)

                # Plot standard deviation bands
                standard_deviation_bands!(cycle, trend, Makie.wong_colors()[j])

                fig.content[1].ylabel = "Mean" 
                if length(gdf) == length(SCENARIOS) 
                    fig.content[i].xlabel = i in [1, 4, 5] ? "Steps" : ""
                    fig.content[i].xticksvisible = i in [2, 3] ? false : true  
                    fig.content[i].xticklabelsvisible = i in [2, 3] ? false : true  
                else
                    fig.content[i].xlabel = "Steps"
                end
                # Add dotted lines to the plot
                add_lines!(gdf[i])
            end
        
            # Set x-axis ticks
            ax.xticks = 0:200:800
        end
        configure_axes_links(fig, gdf)
        fig[end+1,1:3] = Legend(
            fig,
            fig.content[1];
            tellheight = true,
            tellwidth = false,
            orientation = :horizontal
        )
        return fig
    end

    p = plot_figure_5(filter(:shock => ==("Missing"), mdf), [:ion, :iterm], ["ON rate", "Term rate"])
    save("Fig5.pdf", p)

    # Figure 6
    ## Figure 6a
    function plot_figure_6a_C1(df, vars, ylabels)
        fig = Figure(size = (800, 400))
        for i in 1:length(vars)
            ax = fig[1, i] = Axis(fig, title = ylabels[i])
            # Take gdf for dataframe with and without shock
            gdf, gdf_shock = generate_dataframes(df, "Expansionary")
                
            for j in 1:length(gdf)
                # Apply HP filter to the selected variable of the non-shocked df
                _, trend_missing = hp_filter(gdf[j][!, vars[i]], 1600)
                # Apply HP filter to the selected variable of the shocked df
                _, trend_shock = hp_filter(gdf_shock[j][!, vars[i]], 1600)
                # Compute relative impact
                growth_rate_trend = (trend_shock .- trend_missing) ./ trend_missing
                # Plot the main trend line
                lines!(growth_rate_trend, label = only(unique(gdf[j].scenario)), color = j == 1 ? :black : Makie.wong_colors()[j])
            end
            # Set ticks x-axys
            ax.xticks = 0:200:800
            fig.content[1].ylabel = "Mean"
            fig.content[i].xlabel = "Steps"
            linkyaxes!(fig.content...)
            invisible_yaxis!(fig, i)
        end
        fig[end+1,1:length(vars)] = Legend(
            fig,
            fig.content[1];
            tellheight = true,
            tellwidth = false,
            orientation = :horizontal
        )
        return fig
    end

    df1 = generate_df(df, mdf; shock = "Expansionary", banks = true, relative = true);
    p = plot_figure_6a_C1(df1, [:ON_interbank, :Term_interbank], ["Overnight volumes", "Term volumes"])
    save("Fig6a.pdf", p)

    ## Figure 6b
    function plot_figure_6b(dataframes, vars, ylabels)
        fig = Figure(size = (800, 400))
        for i in 1:length(dataframes)
            ax = fig[1, i] = Axis(fig, title = ylabels[i])
            # Take gdf for dataframe with and without shock
            gdf, gdf_shock = generate_dataframes(dataframes, i, "Expansionary")
                
            for j in 1:length(gdf)
                # Apply HP filter to the selected variable of the non-shocked df
                _, trend_missing = hp_filter(gdf[j][!, vars[i]], 1600)
                # Apply HP filter to the selected variable of the shocked df
                _, trend_shock = hp_filter(gdf_shock[j][!, vars[i]], 1600)
                # Compute relative impact
                growth_rate_trend = (trend_shock .- trend_missing) ./ trend_missing
                # Plot the main trend line
                lines!(growth_rate_trend, label = only(unique(gdf[j].scenario)), color = j == 1 ? :black : Makie.wong_colors()[j])
            end
            # Set ticks x-axys
            ax.xticks = 0:200:800
            fig.content[1].ylabel = "Mean"
            fig.content[i].xlabel = "Steps"
            linkyaxes!(fig.content...)
            invisible_yaxis!(fig, i)
        end
        fig[end+1,1:length(dataframes)] = Legend(
            fig,
            fig.content[1];
            tellheight = true,
            tellwidth = false,
            orientation = :horizontal
        )
        return fig
    end

    df_firms = generate_df(df, mdf; shock = "Corridor", relative = true);
    df_hh = generate_df(df, mdf; shock = "Corridor", households = true, relative = true);
    p = plot_figure_6b([df_firms, df_hh], [:GDP, :consumption], ["GDP", "Consumption"])
    save("Fig6b.pdf", p)

    # Figure 7
    function plots_by_group_rationing(df, vars_nominator, vars_denominator)
        fig = Figure(size = (1200, 600))
        gdf = groupby(df, :scenario)

        for j in 1:length(gdf)
            axes = ((1:2,1), (1,2), (1,3), (2,2), (2,3))
            ax = fig[axes[j]...] = Axis(fig, title = SCENARIOS[j])

            for i in eachindex(VALUE_GROUPS)
                # Filter the data for the current status
                sdf = filter(r -> r.value == VALUE_GROUPS[i], gdf[j])

                # Apply HP filter to the selected variable
                cycle, trend = hp_filter((1 .- (sdf[!, vars_nominator[1]]./sdf[!, vars_denominator[1]])) .* 100, 1600)

                # Plot the main trend line
                lines!(trend; linewidth = 2, label = VALUE_GROUPS[i], color = Makie.wong_colors()[2:end][i])

                # Add standard deviation lines
                standard_deviation_bands!(cycle, trend, Makie.wong_colors()[2:end][i])

                # Set ticks x-axys
                ax.xticks = 0:200:800
                fig.content[1].ylabel = "Rate (%)"
                fig.content[j].xlabel = j in [1, 4, 5] ? "Steps" : ""
                fig.content[j].xticksvisible = j in [2, 3] ? false : true  
                fig.content[j].xticklabelsvisible = j in [2, 3] ? false : true  
            end
        end
    configure_axes_links(fig, gdf)
    fig[end+1,1:3] = Legend(
            fig,
            fig.content[1];
            tellheight = true,
            tellwidth = false,
            orientation = :horizontal
        )
        return fig
    end

    df2 = generate_df(df, mdf; shock = "Corridor", banks = true, value = true, status = true);
    p = plots_by_group_rationing(filter(:status => ==("borrower"), df2), [:Term_interbank], [:ib_term_demand])
    save("Fig7.pdf", p)

    #####################################################
    # Appendix
    # Figure C1
    df1 = generate_df(df, mdf; shock = "Corridor-unique", banks = true, relative = true);
    p = plot_figure_6a_C1(df1, [:ON_interbank, :Term_interbank], ["Overnight volumes", "Term volumes"])
    save("FigC1.pdf", p)

    # Figure C2
    p = plots_by_group(df_banks_1, [:Term_interbank])
    save("FigC2.pdf", p)

    # Figure C3 
    p = plots_by_group_rationing(filter(:status => ==("borrower"), df2), [:ON_interbank], [:ib_on_demand])
    save("FigC3.pdf", p)

    # Figure C4
    df_hh = generate_df(df, mdf; shock = "Missing", households = true);
    df_firms = generate_df(df, mdf; shock = "Missing");
    # DEBT to GDP evolution and Public debt
    df_gov = @pipe df |> dropmissing(_, :bills) |> 
        filter(:id => x -> x == 1141, _) |> 
        groupby(_, [:step, :scenario]) |> 
        combine(_, :bills .=> mean, renamecols = false);
    # Add debt_to_gdp variable and compute per scenario per shock => see plots/utils.jl
    debt_to_gdp(df_gov, df_firms)
    p = plot_figures_3c_3d_C4_C5([df_gov], [:debt_to_gdp], ["Debt-to-GDP ratio"])
    save("FigC4.pdf", p)

    # Figure C5
    df_hh = generate_df(df, mdf; shock = "Corridor", households = true);
    df_firms = generate_df(df, mdf; shock = "Corridor");
    # Add debt_to_gdp variable and compute per scenario per shock => see plots/utils.jl
    debt_to_gdp(df_gov, df_firms)
    p = plot_figures_3c_3d_C4_C5([df_gov], [:debt_to_gdp], ["Debt-to-GDP ratio"])
    save("FigC5.pdf", p)

    # Figure C6
    df_firms = generate_df(df, mdf; shock = "Missing");
    p = plot_figure_3a_C6(df_firms, [:GDP, :consumption], ["GDP", "Consumption"])
    save("FigC6.pdf", p)

    # Figure D1
    ## Figure D1a
    df_firms = generate_df(df, mdf; shock = "Expansionary", relative = true);
    df_hh = generate_df(df, mdf; shock = "Expansionary", households = true, relative = true);
    p = plot_figure_6b([df_firms, df_hh], [:GDP, :consumption], ["GDP", "Consumption"])
    save("FigD1a.pdf", p)
    ## Figure D1b
    df_firms = generate_df(df, mdf; shock = "Corridor-unique", relative = true);
    df_hh = generate_df(df, mdf; shock = "Corridor-unique", households = true, relative = true);
    p = plot_figure_6b([df_firms, df_hh], [:GDP, :consumption], ["GDP", "Consumption"])
    save("FigD1b.pdf", p)
end