Greenovation across Boston

Has COVID19 also affected people’s desire to go green?

library(tidyverse)
library(data.table)
library(viridis)
library(hrbrthemes)
library(plotly)
library(gganimate)
library(stringr)
library(lubridate)
library(wrapr)
require(scales)
library(directlabels)
library(tidytext)
library(tm)
library(wordcloud)
library(RColorBrewer)
library(SnowballC)
library(leaflet)
library(sf)
library(tigris)
library(gridExtra)


permit2020 <- fread("D:\\Google Drive\\BARI Research Team Data Library\\Permits\\Outputs\\Permits.Records.Geocoded.SSH.01072020.csv") # 460933
# permit2019 <- fread("D:\\Google Drive\\BARI Research Team Data Library\\Permits\\Data\\Permits.Records.Geocoded.20190805.csv") # 419700
# permit2019 <- permit2019[,c("PermitNumber", "description", "ISSUED_DATE",
#                             "Property_ID","parcel_num", "X" ,"Y" , "Land_Parcel_ID" ,
#                             "Blk_ID_10" , "BG_ID_10" , "CT_ID_10" , "BRA_PD",
#                             "addition", "NSA_Name" )]
# raw2020 <- fread("D:\\Google Drive\\BARI Research Team Data Library\\Permits\\Data\\rawpermits.csv")
# raw2020_1 <- fread("D:/Google Drive/BARI Research Team Data Library/Permits/Outputs/buildingpermitsgeo2020_SSH.csv")

Boston_TA <- fread("C:\\Users\\ux305\\Downloads\\PADCross.Record.2019.csv")
# Boston_TA$CT_ID_10 <- as.character(Boston_TA$CT_ID_10)
# Boston_TA[Boston_TA$CT_ID_10 == "         NA",]$CT_ID_10 <- 
#   landParcels$CT_ID_10[match(Boston_TA[Boston_TA$CT_ID_10 == "         NA",]$Land_Parcel_ID, landParcels$Land_Parcel_ID)]

Building permit dataset opens a window to the citizens’ desires and choices through time. Annually, thousands of citizens decide to renovate their properties and add more energy efficient features to it. Using the building permit of City of Boston, w can have an insight into these trends. Let’s see what categories of green permits exist?
According to the permits description, there are four recognizable categories, and some of them might seem pretty close, but are actually quite different in definition. The categories are: * weatherization
* insulation
* air sealing
* solar cells
Just like any other element tied to human society, buidling permits are also exposed to the pandemic going around. Using the dataset, we can dive in deeper and see how are they distributed across the year and have they all been affected by the pandemic in the same way, if they have been at all?

# add new metrics
permit2020$NOTES <- tolower(permit2020$NOTES)

permit2020$solar <- ifelse(str_detect(permit2020$NOTES, "solar|pv "), 1, 0)
permit2020$insulation <- ifelse(str_detect(permit2020$NOTES, "insul"), 1, 0)
permit2020$weatherization <- ifelse(str_detect(permit2020$NOTES, "weather"), 1, 0)
permit2020$air_sealing <- ifelse(str_detect(permit2020$NOTES, "seal"), 1, 0)

permit2020$greenovation <- ifelse(str_detect(permit2020$NOTES, "insul|solar|pv|energy"), 1, 0)

# subsetting to required columns
permit2020 <- permit2020[,c("PermitNumber", "description", "NOTES", "ISSUED_DATE",
                            "Property_ID","parcel_num", "X" ,"Y" , "Land_Parcel_ID" ,
                            "Blk_ID_10" , "BG_ID_10" , "CT_ID_10" , "BRA_PD", "reno",
                            "addition", "NSA_NAME", "DECLARED_VALUATION", "greenovation",
                            "solar", "insulation", "weatherization", "air_sealing")]
# write.csv()

# filtering for specific desired time period
permit2020$date <- as.Date(permit2020$ISSUED_DATE, "%Y-%m-%d %H:%M:%S", tz = "GMT")
permit2020$month <- month(permit2020$date)
permit2020$Year <- year(permit2020$date)

permit_1920 <- permit2020 %>%
  filter(Year %in% c("2020", "2019")) %>%
  filter(month %in% c(1,2,3,4,5,6))

permit_1920 <- permit_1920 %>%
  filter(reno == 1 | addition == 1)

permit_1920 <- permit_1920 %>%
  group_by(Year, month) %>%
  mutate(freq_permits = n())

permit_1920$month <- factor(month.abb[permit_1920$month],levels=month.abb)

# data ready to use for analysis



permit_1920_green <- permit_1920 %>%
  filter(greenovation == 1)

permit_1920_green <- permit_1920_green %>%
  group_by(Year, month) %>%
  mutate(solar_month_year = sum(solar), insulation_month_year = sum(insulation), 
         weatherization_month_year = sum(weatherization), air_sealing_month_year = sum(air_sealing)) 

p1 <- ggplot(permit_1920_green, aes(month, solar_month_year, fill = as.factor(Year))) + 
  geom_bar(stat = "identity", position=position_dodge()) +
  xlab("") + ylab("Solar cell permits") + # Set axis labels
  ggtitle(label = "Figure 1" ,  subtitle = "Solar cell permits in the first six months of the year") +
  theme_bw() +
  theme(plot.title = element_text(size=8), plot.subtitle = element_text(size = 9)) +
  scale_fill_manual(name  ="Year", values=c("#999999", "#E69F00"))


p2 <- ggplot(permit_1920_green, aes(month, insulation_month_year, fill = as.factor(Year))) + 
  geom_bar(stat = "identity", position=position_dodge()) +
  xlab("") + ylab("Insulation permits") + # Set axis labels
  ggtitle(label = "Figure 2" ,  subtitle = "Insulation permits in the first six months of the year") +
  theme_bw() +
  theme(plot.title = element_text(size=8), plot.subtitle = element_text(size = 9)) +
  scale_fill_manual(name  ="Year", values=c("#999999", "#E69F00"))

p3 <- ggplot(permit_1920_green, aes(month, weatherization_month_year, fill = as.factor(Year))) + 
  geom_bar(stat = "identity", position=position_dodge()) +
  xlab("") + ylab("Weatherization permits") + # Set axis labels
  ggtitle(label = "Figure 3" ,  subtitle = "Weatherization permits in the first six months of the year") +
  theme_bw() +
  theme(plot.title = element_text(size=8), plot.subtitle = element_text(size = 9)) +
  scale_fill_manual(name  ="Year", values=c("#999999", "#E69F00"))

p4 <- ggplot(permit_1920_green, aes(month, air_sealing_month_year, fill = as.factor(Year))) + 
  geom_bar(stat = "identity", position=position_dodge()) +
  xlab("") + ylab("Air sealing permits") + # Set axis labels
  ggtitle(label = "Figure 4" ,  subtitle = "Air sealing permits in the first six months of the year") +
  theme_bw() +
  theme(plot.title = element_text(size=8), plot.subtitle = element_text(size = 9)) +
  scale_fill_manual(name  ="Year", values=c("#999999", "#E69F00"))

grid.arrange(p1, p2, p3, p4)


In this graph, we can see the different assignment of the four gategories of green permits across the first six months of the year. What get the eye is the huge jump in the number of solar cell permits in 2020 before the pandemic. The other three categories also show a head start before the pandemic compared to last year.

In order to test our hypothesis, we need to compare the general trend of issuing green permits from Jan to July in both years of 2019 and 2020. Let’s get into it!
First, let’s see if the trend of issuing all types of permits has changed pre-covid vs during.

# subsetting to required columns
# permit2020 <- permit2020[,c("PermitNumber", "description", "NOTES", "ISSUED_DATE",
#                             "Property_ID","parcel_num", "X" ,"Y" , "Land_Parcel_ID" ,
#                             "Blk_ID_10" , "BG_ID_10" , "CT_ID_10" , "BRA_PD", "reno",
#                             "addition", "NSA_NAME", "DECLARED_VALUATION" )]
# 
# # filtering for specific desired time period
# permit2020$date <- as.Date(permit2020$ISSUED_DATE, "%Y-%m-%d %H:%M:%S", tz = "GMT")
# permit2020$month <- month(permit2020$date)
# permit2020$Year <- year(permit2020$date)
# 
# permit_1920 <- permit2020 %>%
#   filter(Year %in% c("2020", "2019")) %>%
#   filter(month %in% c(1,2,3,4,5,6))

permit_1920 <- permit_1920 %>%
  group_by(Year, month) %>%
  mutate(freq_permits = n())

# permit_1920$month <- factor(month.abb[permit_1920$month],levels=month.abb)

ggplot(data=permit_1920, aes(x=as.factor(month),group=as.factor(Year))) +
  geom_line(stat = "count", aes(color = as.factor(Year)),size=1) +
  scale_color_viridis(discrete = TRUE) +
  ggtitle(label = "Figure 5" ,  subtitle = "Number of permits issued in the first 6 months of the year - 2019 vs 2020") +
  theme_ipsum() +
  ylab("") + xlab("") +
  scale_colour_discrete(name  ="Year") +
  theme(plot.title = element_text(size=8), plot.subtitle = element_text(size = 9))


As apparent in figure 5, the number of permits has droped drastically from Jan to Jun 2020 compared to its trend in 2019, but it still shows an increase starting from Apr. In the month of Apr alone, there is a 77% drop in the total number of permits. This is the highest drop reported yet. 

It would be interesting to look at the growth rate of issuing permits for these two time periods.

# ratio of 2020 to 2019

# permit_1920

dt <- dcast(permit_1920, month ~ Year, value.var = "freq_permits")
names(dt)[2] <- "Year_2019"
names(dt)[3] <- "Year_2020"
dt <- dt %>%
  mutate(ratio20TO19 = Year_2020/Year_2019)

ggplot(data=dt, aes(x=as.factor(month), y = ratio20TO19)) +
  geom_point(shape=21, color="black", fill="Red", size=4) + 
  geom_line(aes(group=1),size=1, color="Dark Blue") +
  ggtitle(label = "Figure 6" ,  subtitle = "Growth rate of permits issued in the first 6 months of the year - 2019 vs 2020") +
  theme_ipsum() +
  ylab("") + xlab("") +
  geom_segment(color="black", 
                  aes(
                    xend=c(tail(month, n=-1), NA), 
                    yend=c(tail(ratio20TO19, n=-1), NA)
                  ),
                  arrow=arrow(length=unit(0.3,"cm"), type = "closed")
      ) +
  theme(plot.title = element_text(size=8), plot.subtitle = element_text(size = 9))


Despite the apparent consistant growth rate in Jan, there is a very visible drop in the following three months. After that, we can see the recovery process begining again with a more steep slope from May to Jun. 

Buidling permits are also an indicator of the residents’ investments in their properties. Let’s see how is that been going during the pandemic in comparison with 2019 when there was no crisis on hand.

permit_1920_1 <- permit_1920 %>%
  group_by(Year, month) %>%
  mutate(Value_tot = sum(DECLARED_VALUATION))

ggplot(data=permit_1920_1, aes(x=as.factor(month), y = Value_tot, group=as.factor(Year))) +
  geom_line(stat = "identity", aes(color = as.factor(Year)),size=1) +
  scale_color_viridis(discrete = TRUE) +
  ggtitle(label = "Figure 7" ,  subtitle = "Total declared value of permits issued in the first 6 months of the year - 2019 vs 2020") +
  theme_ipsum() +
  ylab("") + xlab("") +
  scale_colour_discrete(name  ="Year") +
  theme(plot.title = element_text(size=8), plot.subtitle = element_text(size = 9)) +
  scale_y_continuous(labels = scales::comma)


According to this figure, we also experienced a fall in investment from Feb to March last year. The fall is more drastic this year. Both alteration periods (the Feb to Apr drop and the Apr to Jun rise) have a more steep slope.

Maybe we hould look at the investment rate of 2019 vs 2020 to understand it better.

permit_1920_2 <- permit_1920_1[,c("Year", "month", "Value_tot")]
permit_1920_2 <- permit_1920_2[!duplicated(permit_1920_2),]

dt <- dcast(permit_1920_2, month ~ Year, value.var = "Value_tot")
names(dt)[2] <- "Year_2019"
names(dt)[3] <- "Year_2020"
dt <- dt %>%
  mutate(ratio20TO19 = Year_2020/Year_2019)

ggplot(data=dt, aes(x=as.factor(month), y = ratio20TO19)) +
  geom_point(shape=21, color="black", fill="Red", size=4) + 
  geom_line(aes(group=1),size=1, color="Dark Blue") +
  ggtitle(label = "Figure 8" ,  subtitle = "Growth rate of values of permits issued in the first 6 months of the year - 2019 vs 2020") +
  theme_ipsum() +
  ylab("") + xlab("") +
  geom_segment(color="black", 
                  aes(
                    xend=c(tail(month, n=-1), NA), 
                    yend=c(tail(ratio20TO19, n=-1), NA)
                  ),
                  arrow=arrow(length=unit(0.3,"cm"), type = "closed")
      ) +
  theme(plot.title = element_text(size=8), plot.subtitle = element_text(size = 9))


Here we wintess a more unified fall starting from Jan. This can very much be justified with the anxiety of the job market getting worse as the pandemic found its way into our country.


So far, we analyzed all permits including green ones. But how would all these trends look like if only looking at energy related permits?

# permit_1920_orig <- permit_1920

# permit_1920 <- permit_1920 %>%
#   filter(reno == 1 | addition == 1)
# 
# permit_1920 <- permit_1920 %>%
#   group_by(Year, month) %>%
#   mutate(freq_permits = n())
# 
# permit_1920$NOTES <- tolower(permit_1920$NOTES)
# 
# permit_1920$greenovation <- ifelse(str_detect(permit_1920$NOTES, "insul|solar|pv|energy"), 1, 0)

permit_1920_green <- permit_1920 %>%
  filter(greenovation == 1) 

ggplot(data=permit_1920_green, aes(x=as.factor(month),group=as.factor(Year))) +
  geom_line(stat = "count", aes(color = as.factor(Year)),size=1) +
  scale_color_viridis(discrete = TRUE) +
  ggtitle(label = "Figure 9" ,  subtitle = "Number of green permits issued in the first 6 months of the year - 2019 vs 2020") +
  theme_ipsum() +
  ylab("") + xlab("") +
  scale_colour_discrete(name  ="Year") +
  theme(plot.title = element_text(size=8), plot.subtitle = element_text(size = 9))


In figure 9, we see the exact same pattern we witnessed for all other permits in general with no categories excluded. In the month of Apr again, there is a 77% drop in the total number of permits. This is the highest drop reported yet.  Another cool thing to notice here is that, in figure1, when looking at all permits, two lines counter eachother in early Feb as the fall continues while in figure2, when only looking at green permits, this happens in late Feb. The line actually flattens in Feb after a sharp fall in Jan and then shows an even sharper drop at the start of March. Yet again, the going back to normal process seems to start in early Apr for both graphs.

# ratio of 2020 to 2019

# permit_1920

dt <- dcast(permit_1920_green, month ~ Year, value.var = "freq_permits")
names(dt)[2] <- "Year_2019"
names(dt)[3] <- "Year_2020"
dt <- dt %>%
  mutate(ratio20TO19 = Year_2020/Year_2019)

ggplot(data=dt, aes(x=as.factor(month), y = ratio20TO19)) +
  geom_point(shape=21, color="black", fill="Red", size=4) + 
  geom_line(aes(group=1),size=1, color="Dark Blue") +
  ggtitle(label = "Figure 10" ,  subtitle = "Growth rate of green permits issued in the first 6 months of the year - 2019 vs 2020") +
  theme_ipsum() +
  ylab("") + xlab("") +
  geom_segment(color="black", 
                  aes(
                    xend=c(tail(month, n=-1), NA), 
                    yend=c(tail(ratio20TO19, n=-1), NA)
                  ),
                  arrow=arrow(length=unit(0.3,"cm"), type = "closed")
      ) +
  theme(plot.title = element_text(size=8), plot.subtitle = element_text(size = 9))


Despite what we saw in figure 6 and according to figure 10, there is an immediet drop in the number of green permits at the start of the year. This can be sue to the fact that green permits may be concieved as less essential and more expensive compared to other types of renovations.

Let’s see how the green investments are different in 2020 compared to 2019.

permit_1920_1 <- permit_1920_green %>%
  group_by(Year, month) %>%
  mutate(Value_tot = sum(DECLARED_VALUATION))

ggplot(data=permit_1920_1, aes(x=as.factor(month), y = Value_tot, group=as.factor(Year))) +
  geom_line(stat = "identity", aes(color = as.factor(Year)),size=1) +
  scale_color_viridis(discrete = TRUE) +
  ggtitle(label = "Figure 11" ,  subtitle = "Total declared value of green permits issued in the first 6 months of the year - 2019 vs 2020") +
  theme_ipsum() +
  ylab("") + xlab("") +
  scale_colour_discrete(name  ="Year") +
  theme(plot.title = element_text(size=8), plot.subtitle = element_text(size = 9)) +
  scale_y_continuous(labels = scales::comma)


Figure 11 shows a less apparent trend compared to the other figures featuring COVID19 effect on the building permits. There are two rises, from Feb to Mar and From Apr to May. The other different element of this figure is that there is no sign of a recovery state yet.

Would the growth rate of green value tell the same story?

permit_1920_2 <- permit_1920_1[,c("Year", "month", "Value_tot")]
permit_1920_2 <- permit_1920_2[!duplicated(permit_1920_2),]

dt <- dcast(permit_1920_2, month ~ Year, value.var = "Value_tot")
names(dt)[2] <- "Year_2019"
names(dt)[3] <- "Year_2020"
dt <- dt %>%
  mutate(ratio20TO19 = Year_2020/Year_2019)

ggplot(data=dt, aes(x=as.factor(month), y = ratio20TO19)) +
  geom_point(shape=21, color="black", fill="Red", size=4) + 
  geom_line(aes(group=1),size=1, color="Dark Blue") +
  ggtitle(label = "Figure 12" ,  subtitle = "Growth rate of values of permits issued in the first 6 months of the year - 2019 vs 2020") +
  theme_ipsum() +
  ylab("") + xlab("") +
  geom_segment(color="black", 
                  aes(
                    xend=c(tail(month, n=-1), NA), 
                    yend=c(tail(ratio20TO19, n=-1), NA)
                  ),
                  arrow=arrow(length=unit(0.3,"cm"), type = "closed")
      ) +
  theme(plot.title = element_text(size=8), plot.subtitle = element_text(size = 9))


Figure 12 also presents the same type of pattern where there are two rising points and three drops. More analysis is required to get to the possible reasons for this of course.

The counts are not always the best tools to describe a phenomenon. So, let’s create a metric: what is the percentage of green permits for each parcel?

## calculating the total number of parcels per CT
Boston_TA <- Boston_TA %>%
  group_by(CT_ID_10) %>%
  mutate(num_parcel = n())

permit_1920$CT_ID_10 <- as.character(permit_1920$CT_ID_10)
# permit_1920$month <- factor(month.abb[permit_1920$month],levels=month.abb)


# test 25025000100
# ct = 25025140300 0.12414649
# ct = 25025140400 0.33236394

permit_1920_1 <- permit_1920 %>%
  filter(greenovation == 1) %>%
  filter(!is.na(Land_Parcel_ID)) %>%
  group_by(month, Year, Land_Parcel_ID, CT_ID_10) %>%
  tally() %>%
  mutate(green_parcels = length(unique(Land_Parcel_ID)))

permit_1920_1$num_parcel <- Boston_TA$num_parcel[match(permit_1920_1$CT_ID_10, Boston_TA$CT_ID_10)]
permit_1920_1$metric1 <- (permit_1920_1$green_parcels/permit_1920_1$num_parcel)*100
permit_1920_1 <- permit_1920_1 %>%
  group_by(month, Year) %>%
  mutate(metric_mean = mean(metric1))

ggplot(data= permit_1920_1, aes(x=as.factor(month),y = metric_mean ,group=as.factor(Year))) +
  geom_line(stat = "identity", aes(color = as.factor(Year)),size=1) +
  scale_color_viridis(discrete = TRUE) +
  ggtitle(label = "Figure 13" ,  subtitle = "Percentage of parcels with green permits (per CT) issued in the first 6 months of the year - 2019 vs 2020") +
  theme_ipsum() +
  ylab("") + xlab("") +
  scale_colour_discrete(name  ="Year") +
  theme(plot.title = element_text(size=8), plot.subtitle = element_text(size = 9))


According to figure 13, there is a huge spike in the percentage of parcels with green permits in Mar. Other than this observation, the same pattern is actually visible in 2019; although with a far less aggresive rise. Tihs could be due to a specific project in a census tract. More analysis is needed.

Here’s the fun part! Let’s map the green permits by our new metric and see if any particular part of the city has shown more interest in going green compared to last year.

permit_1920_2 <- permit_1920 %>%
  filter(Year == 2019) %>%
  filter(greenovation == 1) %>%
  filter(!is.na(Land_Parcel_ID)) %>%
  group_by(CT_ID_10) %>%
  mutate(green_parcels = length(unique(Land_Parcel_ID)))

permit_1920_2$num_parcel <- Boston_TA$num_parcel[match(permit_1920_2$CT_ID_10, Boston_TA$CT_ID_10)]
permit_1920_2$metric2 <- (permit_1920_2$green_parcels/permit_1920_2$num_parcel)*100

# mapping

# lookup_code(state="MA",county="Middlesex") # "The code for Massachusetts is '25' and the code for Middlesex County is '017'."
# lookup_code(state="MA",county="Suffolk") # "The code for Massachusetts is '25' and the code for Suffolk County is '025'."
# 
# shapefiles <- tracts(state='25', county=c('17', '25'))
# # plot(shapefiles)

shapefiles <- st_read("C:\\Users\\ux305\\Downloads\\BostonTracts2010-20200630T235634Z-001\\BostonTracts2010\\BostonTracts2010.shp")
## Reading layer `BostonTracts2010' from data source `C:\Users\ux305\Downloads\BostonTracts2010-20200630T235634Z-001\BostonTracts2010\BostonTracts2010.shp' using driver `ESRI Shapefile'
## Simple feature collection with 178 features and 16 fields
## geometry type:  POLYGON
## dimension:      XY
## bbox:           xmin: -71.19115 ymin: 42.22788 xmax: -70.98471 ymax: 42.40493
## CRS:            4326
# mydatamerged <- geo_join(shapefiles, permit_1920_green, "Ln_P_ID", "Land_Parcel_ID")
mydatamerged <- geo_join(shapefiles, permit_1920_2, "CT_ID_10", "CT_ID_10")
mydatamerged$nw_cnst <- ifelse(is.na(mydatamerged$metric2), 0, mydatamerged$metric2)

df <- mydatamerged
# df <- df[!is.na(df$CITY),]
# df <- df[!is.na(df$prop_TR),]
mypopup <- paste0("Parcels with Green Permits: ", ifelse(is.na(df$metric2), 0, round((df$metric2), 2)), "%")
# mypopup <- paste0(df$ISD_Nbhd, ":", "Neighborhood: ", df$, "<br>", "Parcels with Green Permits: ", round((df$metric_mean)*100, 2), "%")
mypal <- colorNumeric(
  palette = "YlGnBu",
  domain = df$metric2
)

# YlOrRd, RdYlGn, YlGnBu, RdYlBu

myLAT <- 42.3398
myLNG <- -71.0892
mycentername <- "NEU"


mymap <- leaflet() %>%
  addProviderTiles("CartoDB.Positron") %>%
  setView(myLNG, myLAT, zoom = 12) %>%
  addMarkers(lat=myLAT, lng=myLNG, popup=mycentername) %>%
  addPolygons(data = df, highlight = highlightOptions(weight = 3,
                                           color = "red",
                                           bringToFront = TRUE) ,
              fillColor = ~mypal(metric2), 
              color = "#000000", # you need to use hex colors - #b2aeae
              fillOpacity = 0.7, 
              weight = 1,
              smoothFactor = 0.2  ,
              popup = mypopup) %>% 

  addLegend(pal = mypal,
            values = df$metric2,
            position = "bottomright",
            title = "Green Permits",
            labFormat = labelFormat(suffix="%", transform = function(x) 1*x))
# %>%
#    addLayersControl(
#     baseGroups = c("OSM (default)", "Toner", "Toner Lite"),
#     overlayGroups = c("Quakes", "Outline"),
#     options = layersControlOptions(collapsed = FALSE)
#   )

mymap


And for 2020:


permit_1920_2 <- permit_1920 %>%
  filter(Year == 2020) %>%
  filter(greenovation == 1) %>%
  filter(!is.na(Land_Parcel_ID)) %>%
  group_by(CT_ID_10) %>%
  mutate(green_parcels = length(unique(Land_Parcel_ID)))

permit_1920_2$num_parcel <- Boston_TA$num_parcel[match(permit_1920_2$CT_ID_10, Boston_TA$CT_ID_10)]
permit_1920_2$metric2 <- (permit_1920_2$green_parcels/permit_1920_2$num_parcel)*100


# mydatamerged <- geo_join(shapefiles, permit_1920_green, "Ln_P_ID", "Land_Parcel_ID")
mydatamerged <- geo_join(shapefiles, permit_1920_2, "CT_ID_10", "CT_ID_10")
mydatamerged$nw_cnst <- ifelse(is.na(mydatamerged$metric2), 0, mydatamerged$metric2)

df <- mydatamerged
# df <- df[!is.na(df$CITY),]
# df <- df[!is.na(df$prop_TR),]
mypopup <- paste0("Parcels with Green Permits: ", ifelse(is.na(df$metric2), 0, round((df$metric2), 2)), "%")
# mypopup <- paste0(df$ISD_Nbhd, ":", "Neighborhood: ", df$, "<br>", "Parcels with Green Permits: ", round((df$metric_mean)*100, 2), "%")
mypal <- colorNumeric(
  palette = "YlGnBu",
  domain = df$metric2
)

# YlOrRd, RdYlGn, YlGnBu, RdYlBu

myLAT <- 42.3398
myLNG <- -71.0892
mycentername <- "NEU"


mymap <- leaflet() %>%
  addProviderTiles("CartoDB.Positron") %>%
  setView(myLNG, myLAT, zoom = 12) %>%
  addMarkers(lat=myLAT, lng=myLNG, popup=mycentername) %>%
  addPolygons(data = df, highlight = highlightOptions(weight = 3,
                                           color = "red",
                                           bringToFront = TRUE) ,
              fillColor = ~mypal(metric2), 
              color = "#000000", # you need to use hex colors - #b2aeae
              fillOpacity = 0.7, 
              weight = 1,
              smoothFactor = 0.2 ,
             popup = mypopup) %>% 
  
  addLegend(pal = mypal,
            values = df$metric2,
            position = "bottomright",
            title = "Green Permits",
            labFormat = labelFormat(suffix="%", transform = function(x) 1*x))
# %>%
#    addLayersControl(
#     baseGroups = c("OSM (default)", "Toner", "Toner Lite"),
#     overlayGroups = c("Quakes", "Outline"),
#     options = layersControlOptions(collapsed = FALSE)
#   )

mymap


There are obvious and also subtle differences between the two maps. The 2019 map is clearly bluer, for example in Dorchester, Hyde Park and East Boston. On the other hand, some tracts that did not have any green permits in the first half of 2019 are now showing progress in going green. In general, unfortunately, there are more tracts with no green permits in 2020. This can have different reasons including getting to a saturation state considering the energy efficiency renovations (you only have to install solar cells once), and of course COVID19. It also might seem like less advantaged neighborhoods have given up on green efforts more that the rest of the city compared to the same time period last year.