I'm trying to iterate through a shapefile with all of the roads in the State of Maine, select only the roads in towns relevant to my project, and then create one new feature class with just those roads. I'm working in Python in a Jupyter notebook and I'm a beginner in both Python and ArcGIS. This code runs, but selects all of the roads, not just those in the town list.
# Create town list
town_list = ['Bangor','Hampden','Veazie','Orono','Brewer','Old Town']
arcpy.env.overwriteOutput = True
for town in town_list:
select = """TOWN = '{}'""".format(town)
arcpy.SelectLayerByAttribute_management(ME_roads,'NEW_SELECTION',select)
arcpy.CopyFeatures_management(ME_roads,relevant_roads)