So, I thought about this a bit. Your problem not a conventional one, and thus it requires a non-obvious solution. The problem is, it's impossible to tailor a solution specifically for it without a good knowledge of your data. So, instead of giving you a set solution, I'll try pointing you to a few, but you'll have to try them yourself and see what works.
The problem. Similarity between lines is not easily definable, but given the configuration of an urban network, I could think of three parameters that affect it:
Azimuth: to avoid one line from being matched to streets perpendicular to it, it's important to filter by the azimuths of the lines. The azimuths can't be a perfect match, because of vectorization (or measuring) imperfections, so a given margin θ is necessary. This margin, however, must not exceed the smallest angle between two connected road segments in your network (if it's a perfect grid, that'd be 90 degrees).
Distance: to avoid two lines with similar azimuths, but in opposite sides of the city, from matching, distance must be taken into account as well. Since traditional distance algorithms only take into account the two closest points between two features, to avoid two connecting streets from matching to each other, I suggest matching both initial and final vertices. This, however, should have a smaller weight than azimuth, because, given that the two shapefiles have lines of different lengths, it's possible that smaller distance be matched to adjascent streets.
- Direction: if you have information on both lanes, this must be taken into account when matching distance and azimuths. If not, you can calculate matching azimuths and counter-azimuths for the same line.
Doing this in QGIS is possible, but quite complicated - it simply wasn't made for this sort of multi-step processing. Since you said you're familiar with python, you can try using GDAL (PyQGIS would also be a possibility, but I'm not as familiar with it).
GDAL is already installed with QGIS. All you'll need to do is download and configure the python bindings for it. To open a shapefile in GDAL:
from osgeo import ogr
driver = ogr.GetDriverByName('ESRI Shapefile')
data_source = driver.Open(<filepath>)
layer = data_source.GetLayer()
for i in range(layer.GetFeatureCount()):
feat = layer.GetFeature(i)
geom = feat.GetGeometryRef()
You can read more about it here. Now, with your data open, you can start doing the checks. Your idea of comparing buffer polygons is not a bad one, and is quite uncomplicated. It may generate incorrect results depending on your network shape, but it's worth giving it a try. You can start by making a search on the neighborhood of your linestring, to limit which features you'll be comparing to. Then checking for intersects:
di_intersect = dict()
buffer_geom1 = geom1.Buffer(radius)
buffer_geom2 = geom2.Buffer(radius)
intersect = buffer_geom1.Intersection(buffer_geom2)
area_intersect = intersect.Area()
feat_id = feat2.GetFieldAsString('FID')
di_intersects[feat_id] = area_intersect
Then you compare which of the features had the largest intersect area, and put that feature's id in a list, or something. Some of these methods might seem confusing if you've never messed with GDAL before, but you can check the documentation here.
If, however, this does not solve the problem, you can try the more laborious way. I can't write the whole code out, but can give you a few snippets to guide your way:
Finding the azimuth (Source)
angle = math.atan2(point2.GetX() - point1.GetX(), point2.GetY() - point1.GetY())
azimuth = math.degrees(angle) if angle > 0 else math.degrees(angle) + 180
Finding closest distance
start1 = geom1.GetPoint(0)
start2 = geom2.GetPoint(0)
end1 = geom1.GetPoint(geom1.GetPointCount() - 1)
end2 = geom2.GetPoint(geom2.GetPointCount() - 1)
dist = start1.Distance(start2) + end1.Distance(end2)
OBS: In the above code, if the features are inverted (i.e. direction doesn't matter), you'll have to match start1 to end2, and start2 to end1.
With this, you can work it in a more complete code, and tweek the search radius to fit your model. The main idea is to iterate through each line in your shapefile 1, and test (either buffer intersection, or azimuth/distance/direction) it for each of the lines in shapefile 2 that are close to it. Then, check which is the best fit, and write that to a list. Final part is to export this list as a CSV, presumibly.
Hope this helps!