python - Python で GeoJSON ポリゴンの面積を取得するにはどうすればよいですか

Question

私はGeoJSONを持っています

{
  "type": "FeatureCollection",
  "features": [
    {
      "type": "Feature",
      "properties": {},
      "geometry": {
        "type": "Polygon",
        "coordinates": [
          [[13.65374516425911, 52.38533382814119], [13.65239769133293, 52.38675829106993], [13.64970274383571, 52.38675829106993], [13.64835527090953, 52.38533382814119], [13.64970274383571, 52.38390931824483], [13.65239769133293, 52.38390931824483], [13.65374516425911, 52.38533382814119]]
        ]
      }
    }
  ]
}

http://geojson.ioは次のように表示されます

その面積(87106.33m^2)をPythonで計算したいと思います。それ、どうやったら出来るの？

私が試したこと

# core modules
from functools import partial

# 3rd pary modules
from shapely.geometry import Polygon
from shapely.ops import transform
import pyproj

l = [[13.65374516425911, 52.38533382814119, 0.0], [13.65239769133293, 52.38675829106993, 0.0], [13.64970274383571, 52.38675829106993, 0.0], [13.64835527090953, 52.38533382814119, 0.0], [13.64970274383571, 52.38390931824483, 0.0], [13.65239769133293, 52.38390931824483, 0.0], [13.65374516425911, 52.38533382814119, 0.0]]
polygon = Polygon(l)

print(polygon.area)
proj = partial(pyproj.transform, pyproj.Proj(init='epsg:4326'),
                   pyproj.Proj(init='epsg:3857'))
print(transform(proj, polygon).area)

そして1.1516745933889345e-05、233827.03300877335最初のものは意味をなさないことが予想されていましたが、2番目のものを修正するにはどうすればよいですか? pyproj.Proj（ initパラメータの設定方法がわかりません）

私epsg:4326はそれがWGS84（ソース）であると推測していますが、epsg:3857については不明です。

より良い結果

以下はもっと近いです：

# core modules
from functools import partial

# 3rd pary modules
import pyproj
from shapely.geometry import Polygon
import shapely.ops as ops

l = [[13.65374516425911, 52.38533382814119, 0],
     [13.65239769133293, 52.38675829106993, 0],
     [13.64970274383571, 52.38675829106993, 0],
     [13.64835527090953, 52.38533382814119, 0],
     [13.64970274383571, 52.38390931824483, 0],
     [13.65239769133293, 52.38390931824483, 0],
     [13.65374516425911, 52.38533382814119, 0]]
polygon = Polygon(l)

print(polygon.area)
geom_area = ops.transform(
    partial(
        pyproj.transform,
        pyproj.Proj(init='EPSG:4326'),
        pyproj.Proj(
            proj='aea',
            lat1=polygon.bounds[1],
            lat2=polygon.bounds[3])),
    polygon)
print(geom_area.area)

それは 87254.7m^2 を与えます - それは geojson.io が言うこととはまだ 148m^2 異なります。なぜそうなのですか？