# ---
# jupyter:
#   jupytext:
#     text_representation:
#       extension: .py
#       format_name: light
#       format_version: '1.5'
#       jupytext_version: 1.14.4
#   kernelspec:
#     display_name: Python 3
#     language: python
#     name: python3
# ---

# # Grab USDA county-level maize data

# +
import requests
import pandas as pd
import numpy as np
import geopandas as gp

api_key = '210BA222-FC6E-3FB2-B4D7-DA2DAA1CC829'

# +
state_alphas = ["AL", "AK", "AZ", "AR", "CA", "CO", "CT", "DE", "FL", "GA",
                "HI", "ID", "IL", "IN", "IA", "KS", "KY", "LA", "ME", "MD",
                "MA", "MI", "MN", "MS", "MO", "MT", "NE", "NV", "NH", "NJ",
                "NM", "NY", "NC", "ND", "OH", "OK", "OR", "PA", "RI", "SC",
                "SD", "TN", "TX", "UT", "VT", "VA", "WA", "WV", "WI", "WY"]

needless_info = ["CV (%)", "agg_level_desc", "asd_desc", "begin_code", "asd_code", 
                 "class_desc", "commodity_desc", "congr_district_code", "state_name",
                 "country_code", "country_name", "domain_desc",
                 "domaincat_desc", "end_code", "freq_desc", "group_desc",
                 "load_time", "location_desc", "prodn_practice_desc",
                 "reference_period_desc", "region_desc", "sector_desc",
                 "short_desc", "source_desc", "statisticcat_desc", "unit_desc",
                 "util_practice_desc", "watershed_code", "watershed_desc",
                 "week_ending", "zip_5", "county_ansi", "state_ansi"]


def get_corn(states, yields=True):
    """
    Grabs county-level corn data from USDA API.
    Input: states = list of U.S. state codes
           yields = boolean, if True return yields in bushels/acre, else return area in acres
    Output: pandas dataframe
    """
    data = pd.DataFrame()
    for state in states:
        print("Now grabbing: " + state)
        if yields:
            dat = requests.get("http://quickstats.nass.usda.gov/api/api_GET/?key=" + api_key + "&source_desc=SURVEY&sector_desc=CROPS&group_desc=FIELD CROPS&commodity_desc=CORN&statisticcat_desc=YIELD&util_practice_desc=GRAIN&unit_desc=BU / ACRE&agg_level_desc=COUNTY&year__GE=1950&state_alpha=" + state)
        else:
            dat = requests.get("http://quickstats.nass.usda.gov/api/api_GET/?key=" + api_key + "&source_desc=SURVEY&sector_desc=CROPS&group_desc=FIELD CROPS&commodity_desc=CORN&statisticcat_desc=AREA HARVESTED&util_practice_desc=GRAIN&unit_desc=ACRES&agg_level_desc=COUNTY&year__GE=1950&state_alpha=" + state)
        if dat.status_code == 200:
            print("Data grabbed from USDA successfully...")
            dat = dat.json()
            dat = pd.DataFrame(dat["data"])
            data = pd.concat([data, dat], ignore_index=True)
            print("Filled!")
        else:
            print("Data grabbed from USDA unsuccessfully. Error code "
                  + str(dat.status_code) + ". Skipping.")
    return data.drop(columns=needless_info)


# + jupyter={"outputs_hidden": true}
# Get data
usda_county_yield = get_corn(state_alphas)
print('\n\n\n Now for the areas...\n\n\n')
usda_county_area = get_corn(state_alphas, yields=False)

# +
# Merge and tidy
usda_county_yield.rename(columns = {'Value':'yield', 'state_fips_code':'state'}, inplace=True)
usda_county_yield['fips'] = usda_county_yield['state'] + usda_county_yield['county_code']
usda_county_yield.drop(columns = ['county_code'], inplace=True)

usda_county_area.rename(columns = {'Value':'area', 'state_fips_code':'state'}, inplace=True)
usda_county_area['fips'] = usda_county_area['state'] + usda_county_area['county_code']
usda_county_area.drop(columns = ['county_code'], inplace=True)

usda_county_all = pd.merge(usda_county_yield.drop_duplicates(subset=['fips','year']).query('county_name != "OTHER (COMBINED) COUNTIES"'),
                           usda_county_area.drop_duplicates(subset=['fips','year']).query('county_name != "OTHER (COMBINED) COUNTIES"'),
                           how='inner')
# -

# County shapefile with lat, lon of centroid
county_shapefile = gp.read_file('../other/plotting_tools/counties_contig_plot.shp')
county_shapefile_lat_lon = county_shapefile.to_crs("EPSG:4326")
county_shapefile_lat_lon['lat'] = county_shapefile_lat_lon['geometry'].geometry.centroid.y
county_shapefile_lat_lon['lon'] = county_shapefile_lat_lon['geometry'].geometry.centroid.x

# Merge
county_all = pd.merge(county_shapefile_lat_lon.drop(columns='geometry'), usda_county_all, on = 'fips', how='inner')

# To numeric
county_all['yield'] = county_all['yield'].astype(float)
county_all['area'] = county_all['area'].apply(lambda x: float(x.replace(",","")) if type(x) == str else x)
# Log yield
county_all['log_yield'] = county_all['yield'].apply(lambda x: np.log(x) if x > 0.0 else np.nan)

county_all

# Save
county_all.to_csv("../data/usda/maize_county_yield_area.csv", index=False)

import pandas as pd
county_all = pd.read_csv("../data/usda/maize_county_yield_area.csv")

county_all = county_all.query('year >= 1956 and year <= 2005')

county_all = county_all.query('lon >= -100')

(county_all.groupby('fips').count()['log_yield'].values == 50).sum()

county_all

# +
############# NOT NEEDED
# # Build index structure for final dataframe (including years)
# all_fips = county_all['fips'].unique()
# years = np.arange(1950, 2021, 1)

# all_fips = [[fips]*len(years) for fips in all_fips]
# all_fips = np.ndarray.flatten(np.asarray(all_fips))

# years = [years] * len(all_fips)
# years = np.ndarray.flatten(np.asarray(years))

# tuples = list(zip(*[all_fips, years]))

# index = pd.MultiIndex.from_tuples(tuples)

# # Build empty dataframe with complete indexing
# all_index = pd.DataFrame(index = index)
# all_index.index.names = ["fips", "year"]

# # Merge yields so that empty county/year pairings are "NaN"
# county_all_indexed = pd.merge(all_index, county_all, on = ['fips', 'year'],  how = "outer")

StarcoderdataJupyterScriptsDedupFiltered