# -*- coding: utf-8 -*-

# ------------------------------------------------------------------------------
# # Yük tabanı mı, denizcilik yetkinliği mi? Ülkelerin denizcilikteki konumunu iki katmanda ölçen bir endeks, sıra belirsizliği ve Türkiye için kaldıraç analizi, 2014–2025
#
# *Cargo base or maritime capability? A two-layer index of countries' maritime position, rank uncertainty and a leverage analysis for Türkiye, 2014–2025*
#
# **GOMD Yayınlar · Deniz ekonomisi 01** — https://gomd.pirireis.edu.tr/publications.php?id=me01
#
# Bu not defteri makaledeki **tüm sayıları, tabloları ve şekilleri** üreten kodun Colab sürümüdür. Veri, replikasyon paketinde (11 Kasım 2026 indirmesi) varsa oradan; yoksa özgün sağlayıcıdan (UNCTADstat) indirilir. Bütün sayılar `out/results.json` dosyasına yazılır.
#
# This notebook reproduces **every number, table and figure** of the paper. Data are taken from the replication package (retrieved 11 November 2026) if present; otherwise they are downloaded from the original provider. All numbers are written to `out/results.json`.
#
# *Makale hazırlık aşamasındadır; dergi adı ve tam metin kabulden sonra eklenecektir. / The paper is in preparation; journal name and full text will be added after acceptance.*
# ------------------------------------------------------------------------------

# ------------------------------------------------------------------------------
# ## 0. Kurulum / Setup
# ------------------------------------------------------------------------------
# [Colab Magic] !pip install -q pandas numpy statsmodels scipy matplotlib adjustText py7zr

# ------------------------------------------------------------------------------
# ## 1. Veri / Data
# Exact reproduction in Colab: upload the replication package (.zip) from the GOMD Publications page, unzip it, then run the cells.
# ------------------------------------------------------------------------------
VERI_KAYNAGI = 'kaynak'   # Colab'da 'zip' secenegi icin not defterine bakin / see the notebook for the 'zip' option in Colab
# Veri / Data. Paketteki veri varsa onu kullanır (11 Kasım 2026 indirmesi, birebir yeniden üretim);
# yoksa özgün kaynaktan indirir (güncel sürüm; sağlayıcı revizyonları nedeniyle sonuçlar küçük farklar gösterebilir).
# Uses the data shipped in the replication package if present (retrieved 11 Nov 2026, exact reproduction);
# otherwise downloads from the original provider (current release; results may differ slightly after provider revisions).
import os, glob, json, urllib.request
import pandas as pd
os.makedirs('data', exist_ok=True); os.makedirs('out', exist_ok=True); os.makedirs('fig', exist_ok=True)
def _unctad(code):
    d = 'data/US_' + code
    if glob.glob(d + '/*.csv'): return
    arc = d + '.7z'
    if not os.path.exists(arc):
        print('UNCTADstat indiriliyor / downloading:', code)
        req = urllib.request.Request('https://unctadstat-api.unctad.org/bulkdownload/US.%s/US_%s' % (code, code), headers={'User-Agent': 'gomd-replication/1.0'})
        with urllib.request.urlopen(req, timeout=300) as r, open(arc, 'wb') as f: f.write(r.read())
    import py7zr
    with py7zr.SevenZipFile(arc) as z: z.extractall(d)
for _c in ['ShipBuilding', 'ShipScrapping', 'FleetBeneficialOwners', 'MerchantFleet', 'SeaborneTrade', 'ContPortThroughput', 'Seafarers', 'VesselValueByOwnership', 'GDPTotal', 'LSCI']: _unctad(_c)

# ------------------------------------------------------------------------------
# ## 2. Göstergeler ve dünya payları / Indicators and world shares
# ------------------------------------------------------------------------------
# Makale 1: veri hazırlığı. UNCTADstat serilerinden ülke x gösterge dünya payları (%), iki dönem ve yıllık
import pandas as pd, numpy as np, glob
D='data/'
def L(n): return pd.read_csv(glob.glob(D+n+'/*.csv')[0],low_memory=False,dtype={'Economy':str,'BeneficialOwnership':str,'FlagOfRegistration':str,'CargoType':str,'ShipType':str,'SeafarerType':str})
def piv(d,code,val,world='0000'):
    e=d[d[code].str.len()==3].pivot_table(index=code,columns='Year',values=val,aggfunc='sum'); w=d[d[code]==world].groupby('Year')[val].sum(); return e,w
sb=L('US_ShipBuilding'); B,Bw=piv(sb,'Economy','Gross Tonnage')
sc=L('US_ShipScrapping'); S,Sw=piv(sc,'Economy','Gross Tonnage')
bo=L('US_FleetBeneficialOwners'); t=bo[bo.FlagOfRegistration=='0000']                      # tüm bayraklar
assert len(t)>0; O,Ow=piv(t,'BeneficialOwnership','Dead weight tons in thousands')
nat=bo[(bo.FlagOfRegistration=='C00001')&(bo.BeneficialOwnership.str.len()==3)].pivot_table(index='BeneficialOwnership',columns='Year',values='Dead weight tons in thousands',aggfunc='sum')   # kendi bayrağında kendi sahipli
mf=L('US_MerchantFleet'); t=mf[mf.ShipType.isin(['0','00'])]; assert len(t)>0; Fl,Fw=piv(t,'Economy','Dead weight tons in thousands')
st=L('US_SeaborneTrade'); T_,Tw=piv(st[st.CargoType=='30'],'Economy','Metric tons in thousands')
cp=L('US_ContPortThroughput'); C_,Cw=piv(cp,'Economy','TEU (Twenty foot Equivalent Unit)')
sf=L('US_Seafarers'); sf.loc[sf.Economy=='579','Economy']='578'; K_,Kw=piv(sf[sf.SeafarerType=='0'],'Economy','Absolute value'); Ko,_=piv(sf[sf.SeafarerType=='1'],'Economy','Absolute value')
vv=L('US_VesselValueByOwnership'); V_=vv[vv.BeneficialOwnership.str.len()==3].pivot_table(index='BeneficialOwnership',columns='Year',values='Percentage of global fleet value')
gd=L('US_GDPTotal'); gc=[c for c in gd.columns if c.startswith('US$ per capita at current')][0]; GP=gd[gd.Economy.str.len()==3].pivot_table(index='Economy',columns='Year',values=gc); GT=gd[gd.Economy.str.len()==3].pivot_table(index='Economy',columns='Year',values='Millions of US$ at current prices')
ls=L('US_LSCI'); ls['Year']=ls.Quarter.str[:4].astype(int); LS=ls[ls.Economy.str.len()==3].pivot_table(index='Economy',columns='Year',values='Index (average Q1 2023 = 100)')
lab={}
for d,c in [(st,'Economy'),(mf,'Economy'),(bo,'BeneficialOwnership')]:
    for k,v in d[[c,c+' Label']].drop_duplicates().values: lab.setdefault(k,v)
lab['792']='Türkiye'
def sh(E,W,yrs,how):                       # dönem payı (%): toplamların oranı ya da yıllık payların ortalaması
    yrs=[y for y in yrs if y in E.columns]
    return (E[yrs].sum(axis=1)/W[yrs].sum()*100) if how=='sum' else (E[yrs].div(W[yrs],axis=1)*100).mean(axis=1)
PER={'E':dict(flow=[2014,2015,2016],stock=[2015,2016,2017],trade=[2014,2015,2016],crew=[2015],val=[2019],gdp=[2014,2015,2016]),
     'L':dict(flow=[2023,2024,2025],stock=[2024,2025,2026],trade=[2022,2023,2024],crew=[2021],val=[2024,2025,2026],gdp=[2022,2023,2024])}
def table(p):
    q=PER[p]; idx=T_.index
    X=pd.DataFrame({'insa':sh(B,Bw,q['flow'],'sum'),'sahiplik':sh(O,Ow,q['stock'],'mean'),'bayrak':sh(Fl,Fw,q['stock'],'mean'),'gemiadami':sh(K_,Kw,q['crew'],'mean'),'geri_donusum':sh(S,Sw,q['flow'],'sum'),
                    'yuk':sh(T_,Tw,q['trade'],'sum'),'teu':sh(C_,Cw,q['trade'],'sum'),'ulusal_bayrak':sh(nat,Fw,q['stock'],'mean'),'filo_degeri':V_[[y for y in q['val'] if y in V_.columns]].mean(axis=1),
                    'zabit':sh(Ko,sf[(sf.Economy=='0000')&(sf.SeafarerType=='1')].groupby('Year')['Absolute value'].sum(),q['crew'],'mean')}).reindex(idx).fillna(0)
    X['gdp_pc']=GP[q['gdp']].mean(axis=1).reindex(idx); X['gdp']=GT[q['gdp']].mean(axis=1).reindex(idx); X['lsci']=LS[[y for y in q['gdp'] if y in LS.columns]].mean(axis=1).reindex(idx); X['ad']=[lab[i] for i in idx]
    return X
CAP=['insa','sahiplik','bayrak','gemiadami']; EXT=CAP+['geri_donusum']; CAR=['yuk','teu']

TRN={"China":"Çin","Bangladesh":"Bangladeş","India":"Hindistan","Republic of Korea":"Kore Cumhuriyeti","Japan":"Japonya","Greece":"Yunanistan","Liberia":"Liberya","Singapore":"Singapur","Panama":"Panama","Marshall Islands":"Marshall Adaları","China, Hong Kong SAR":"Hong Kong (Çin)","Türkiye":"Türkiye","Pakistan":"Pakistan","Denmark":"Danimarka","Malta":"Malta","Indonesia":"Endonezya","Germany":"Almanya","Norway":"Norveç","Italy":"İtalya","Finland":"Finlandiya","France":"Fransa","Russian Federation":"Rusya","Netherlands (Kingdom of the)":"Hollanda","Romania":"Romanya","Poland":"Polonya","Spain":"İspanya","Croatia":"Hırvatistan","Philippines":"Filipinler","United States":"ABD","United Kingdom":"Birleşik Krallık","Ukraine":"Ukrayna","Viet Nam":"Vietnam","Brazil":"Brezilya","China, Taiwan Province of":"Tayvan (Çin)","Bahamas":"Bahamalar","United Arab Emirates":"BAE","Cyprus":"Kıbrıs","Malaysia":"Malezya","Australia":"Avustralya","Saudi Arabia":"Suudi Arabistan","Belgium":"Belçika","Canada":"Kanada","Mexico":"Meksika","Thailand":"Tayland","Switzerland":"İsviçre","Myanmar":"Myanmar","Cuba":"Küba","Oman":"Umman","Morocco":"Fas","South Africa":"Güney Afrika","Colombia":"Kolombiya","Egypt":"Mısır","Iran (Islamic Republic of)":"İran","Portugal":"Portekiz","Sweden":"İsveç","Bermuda":"Bermuda","Monaco":"Monako","Nigeria":"Nijerya","Qatar":"Katar","Kuwait":"Kuveyt","Bulgaria":"Bulgaristan","Sri Lanka":"Sri Lanka","Chile":"Şili","Peru":"Peru","Argentina":"Arjantin","Israel":"İsrail","Algeria":"Cezayir","Iraq":"Irak","Kazakhstan":"Kazakistan","Isle of Man":"Man Adası","Antigua and Barbuda":"Antigua ve Barbuda","Latvia":"Letonya","Estonia":"Estonya","Lithuania":"Litvanya","Ireland":"İrlanda","Georgia":"Gürcistan","Honduras":"Honduras","Ghana":"Gana","Madagascar":"Madagaskar","Cambodia":"Kamboçya","Tuvalu":"Tuvalu","Cameroon":"Kamerun","Saint Vincent and the Grenadines":"Saint Vincent ve Grenadinler","Belize":"Belize","Sierra Leone":"Sierra Leone","Palau":"Palau","Comoros":"Komorlar","Togo":"Togo","Cook Islands":"Cook Adaları","Korea, Dem. People's Rep. of":"Kuzey Kore","Tanzania, United Republic of":"Tanzanya","Venezuela (Bolivarian Rep. of)":"Venezuela","Angola":"Angola","Azerbaijan":"Azerbaycan","New Zealand":"Yeni Zelanda","Dominican Republic":"Dominik Cumhuriyeti","Ecuador":"Ekvador","Kenya":"Kenya","Jamaica":"Jamaika","Guatemala":"Guatemala","Lebanon":"Lübnan","Tunisia":"Tunus","Libya":"Libya","Syrian Arab Republic":"Suriye","Jordan":"Ürdün","Bahrain":"Bahreyn","Kiribati":"Kiribati","Vanuatu":"Vanuatu","Gabon":"Gabon","Barbados":"Barbados","Saint Kitts and Nevis":"Saint Kitts ve Nevis","Mongolia":"Moğolistan","Gibraltar":"Cebelitarık","Cayman Islands":"Cayman Adaları","Luxembourg":"Lüksemburg","Austria":"Avusturya","Czechia":"Çekya","Montenegro":"Karadağ","Albania":"Arnavutluk","Slovenia":"Slovenya","Iceland":"İzlanda","Papua New Guinea":"Papua Yeni Gine","Mozambique":"Mozambik","Guinea":"Gine","Senegal":"Senegal","Côte d'Ivoire":"Fildişi Sahili","Djibouti":"Cibuti","Sudan":"Sudan","Yemen":"Yemen","Uruguay":"Uruguay","Trinidad and Tobago":"Trinidad ve Tobago","Brunei Darussalam":"Brunei","Maldives":"Maldivler","Fiji":"Fiji","Mauritius":"Mauritius","Moldova, Republic of":"Moldova","Turkmenistan":"Türkmenistan","Guinea-Bissau":"Gine-Bissau","Equatorial Guinea":"Ekvator Ginesi","Congo":"Kongo","Dem. Rep. of the Congo":"Kongo DC","Eswatini":"Esvatini","Dominica":"Dominika","Sao Tome and Principe":"Sao Tome ve Principe","Micronesia (Federated States of)":"Mikronezya","Faroe Islands":"Faroe Adaları","Guyana":"Guyana","San Marino":"San Marino","Niue":"Niue","Hungary":"Macaristan","Slovakia":"Slovakya","Paraguay":"Paraguay","Ethiopia":"Etiyopya","Curaçao":"Curaçao","Namibia":"Namibya","Mauritania":"Moritanya","Benin":"Benin","Liberia ":"Liberya"}

# ------------------------------------------------------------------------------
# ## 3. Endeksler, belirsizlik, regresyon, kaldıraç ve şekiller / Indices, uncertainty, regression, leverage and figures
# ------------------------------------------------------------------------------
# Makale 1: çekirdek denizcilik yetkinliği (inşa, sahiplik, bayrak, gemiadamı) ve yük tabanı endeksleri; genişletilmiş (geri dönüşüm dahil) seçenek;
# sıra belirsizliği (Monte Carlo), yoğunluk regresyonu, değişimin gösterge katkıları, kaldıraç analizi (UNCTADstat)
import json, statsmodels.api as sm, matplotlib
matplotlib.use('Agg'); import matplotlib.pyplot as plt
O_='out/'; F_='fig/'; R={}; TR='792'
XE,XL=table('E'),table('L')
for X in (XE,XL):
    X['ad']=[TRN.get(a,a) for a in X.ad]
    X['C']=X[CAP].mean(axis=1); X['C5']=X[EXT].mean(axis=1); X['Y']=X[CAR].mean(axis=1); X['yog']=X.C/X.Y.where(X.Y>0)
    for k in ['C','C5','Y']+EXT+CAR: X['r_'+k]=X[k].rank(ascending=False,method='min')
XL.round(4).to_csv(O_+'endeks_2023_25.csv'); XE.round(4).to_csv(O_+'endeks_2014_16.csv')
R['n_eco']=int(len(XL)); R['n_pos']=int((XL.C>0).sum())
R['sum_check']={k:round(float(XL[k].sum()),2) for k in EXT+CAR}
def row(X,i): r=X.loc[i]; return {k:(r[k] if k=='ad' else (None if pd.isna(r[k]) else round(float(r[k]),4))) for k in X.columns}
R['tr']={'E':row(XE,TR),'L':row(XL,TR)}
top=XL.sort_values('C',ascending=False).head(30).index; R['top30']=[dict(kod=i,**row(XL,i),C_E=round(float(XE.C[i]),4),rC_E=int(XE.r_C[i]),Y_E=round(float(XE.Y[i]),4)) for i in top]
R['hhi']={p:{k:round(float((X[k]**2).sum()),0) for k in EXT+CAR} for p,X in (('E',XE),('L',XL))}
R['top5share']={k:round(float(XL[k].nlargest(5).sum()),1) for k in EXT+CAR}
R['spearman']=XL[XL.C5>0.01][EXT+CAR].corr(method='spearman').round(2).to_dict(); R['spearman_n']=int((XL.C5>0.01).sum())
# 1) Sıra belirsizliği
rng=np.random.default_rng(20261111); N=20000; S_=XL[XL.C5>=0.05].index; R['mc_n']=int(len(S_)); R['mc_N']=N; ranks=np.zeros((N,len(S_))); jt=list(S_).index(TR); opt=np.zeros((N,6),dtype=int)
for n in range(N):
    Z=XL.loc[S_,CAP].copy(); o=(rng.random(6)<.5).astype(int)
    if o[0]: Z['geri_donusum']=XL.loc[S_,'geri_donusum']
    if o[1]: Z['bayrak']=XL.loc[S_,'ulusal_bayrak']
    if o[2]: Z['sahiplik']=XL.loc[S_,'filo_degeri']
    if o[3]: Z['gemiadami']=XL.loc[S_,'zabit']
    m=Z.shape[1]; w=rng.dirichlet(np.ones(m)) if o[4] else np.ones(m)/m
    v=(Z.values*w).sum(1) if o[5] else ((np.sqrt(Z.values)*w).sum(1))**2
    ranks[n]=(-v).argsort().argsort()+1; opt[n]=o
MC=pd.DataFrame({'ad':XL.ad[S_],'C':XL.C[S_],'r0':XL.C[S_].rank(ascending=False,method='min'),'r5':XL.C5[S_].rank(ascending=False,method='min'),'med':np.median(ranks,0),'p05':np.percentile(ranks,5,axis=0),'p95':np.percentile(ranks,95,axis=0)},index=S_).sort_values('C',ascending=False)
MC.round(3).to_csv(O_+'sira_belirsizligi.csv'); R['mc']=[dict(kod=i,ad=r.ad,C=round(r.C,3),r0=int(r.r0),r5=int(r.r5),med=float(r.med),p05=float(r.p05),p95=float(r.p95)) for i,r in MC.head(30).iterrows()]
rt=ranks[:,jt]; R['mc_tr']=dict(r0=int(MC.r0[TR]),r5=int(MC.r5[TR]),med=float(np.median(rt)),p05=float(np.percentile(rt,5)),p95=float(np.percentile(rt,95)),top20=round(float((rt<=20).mean()*100),1),top15=round(float((rt<=15).mean()*100),1),top10=round(float((rt<=10).mean()*100),1))
nm=['geri_donusum_dahil','ulusal_bayrak','filo_degeri','zabit','rastgele_agirlik','aritmetik']
R['mc_tr_by_opt']={nm[k]:dict(evet=float(np.median(rt[opt[:,k]==1])),hayir=float(np.median(rt[opt[:,k]==0]))) for k in range(6)}
# Türkiye sırasındaki varyansın seçeneklere dağılımı (birinci derece duyarlılık: Var(E[r|seçenek])/Var(r))
R['mc_tr_S1']={nm[k]:round(float(np.var([rt[opt[:,k]==a].mean() for a in (0,1)])/rt.var()),3) for k in range(6)}
def rk(Z,gm=False):
    v=Z.mean(axis=1) if not gm else (np.sqrt(Z).mean(axis=1))**2; return int(v.rank(ascending=False,method='min')[TR])
Zb=XL.loc[S_,CAP]; alt={'temel':rk(Zb),'karekok':rk(Zb,True),'genisletilmis':rk(XL.loc[S_,EXT])}
for k in CAP: alt['cikar_'+k]=rk(Zb.drop(columns=k))
for a,b in (('bayrak','ulusal_bayrak'),('sahiplik','filo_degeri'),('gemiadami','zabit')):
    Z=Zb.copy(); Z[a]=XL.loc[S_,b]; alt[b]=rk(Z)
R['alt_rank_tr']=alt
ZbE=XE.loc[S_,CAP]; R['alt_rank_tr_E']={'temel':int(XE.C[S_].rank(ascending=False,method='min')[TR]),'genisletilmis':int(XE.C5[S_].rank(ascending=False,method='min')[TR]),'cikar_gemiadami':int(ZbE.drop(columns='gemiadami').mean(axis=1).rank(ascending=False,method='min')[TR])}
# 2) Yoğunluk regresyonu
def reg(X,col):
    d=X[(X[col]>0)&(X.Y>=0.02)&X.gdp_pc.notna()]; y=np.log(d[col]); Xm=sm.add_constant(pd.DataFrame({'lnY':np.log(d.Y),'lnG':np.log(d.gdp_pc)})); m=sm.OLS(y,Xm).fit(cov_type='HC3'); res=m.resid; rr=res.rank(ascending=False)
    return dict(n=int(len(d)),b=round(float(m.params.lnY),3),b_se=round(float(m.bse.lnY),3),c=round(float(m.params.lnG),3),c_se=round(float(m.bse.lnG),3),c_p=round(float(m.pvalues.lnG),4),r2=round(float(m.rsquared),3),p_b1=round(float(m.t_test('lnY = 1').pvalue),4),tr_res=round(float(res[TR]),3),tr_res_rank=int(rr[TR]),tr_mult=round(float(np.exp(res[TR])),2)),res
R['reg']={}
for p,X in (('E',XE),('L',XL)):
    R['reg'][p],res=reg(X,'C'); X['artik']=res; R['reg'][p+'5'],_=reg(X,'C5')
d=XL[XL.artik.notna()].sort_values('artik',ascending=False); f=lambda q:[dict(ad=r.ad,res=round(float(r.artik),2),C=round(float(r.C),3),Y=round(float(r.Y),3)) for _,r in q.iterrows()]
R['res_top']=f(d.head(12)); R['res_bot']=f(d[d.Y>=0.5].tail(12)); R['res_big']=f(d[d.Y>=1.0])
R['lsci_corr']=dict(C=round(float(XL[XL.lsci.notna()][['C','lsci']].corr(method='spearman').iloc[0,1]),2),Y=round(float(XL[XL.lsci.notna()][['Y','lsci']].corr(method='spearman').iloc[0,1]),2),n=int(XL.lsci.notna().sum()))
# 3) Değişimin gösterge katkıları
def contrib(i): return {k:round(float((XL.loc[i,k]-XE.loc[i,k])/4),4) for k in CAP}
def contrib3(i): return {k:round(float((XL.loc[i,k]-XE.loc[i,k])/3),4) for k in ['insa','sahiplik','bayrak']}
peers=['792','300','380','724','578','208','528','276','826','250','410','356','360','643','076','484','704','818']
pr=lambda i:dict(kod=i,ad=XL.ad[i],C_E=round(float(XE.C[i]),3),C_L=round(float(XL.C[i]),3),C5_E=round(float(XE.C5[i]),3),C5_L=round(float(XL.C5[i]),3),Y_E=round(float(XE.Y[i]),3),Y_L=round(float(XL.Y[i]),3),rC_E=int(XE.r_C[i]),rC_L=int(XL.r_C[i]),rC5_L=int(XL.r_C5[i]),rY_E=int(XE.r_Y[i]),rY_L=int(XL.r_Y[i]),yog_E=round(float(XE.yog[i]),3),yog_L=round(float(XL.yog[i]),3),art=None if pd.isna(XL.artik.get(i,np.nan)) else round(float(XL.artik[i]),2),kat=contrib(i),kat3=contrib3(i),C3_E=round(float(XE.loc[i,['insa','sahiplik','bayrak']].mean()),3),C3_L=round(float(XL.loc[i,['insa','sahiplik','bayrak']].mean()),3),**{k:round(float(XL.loc[i,k]),3) for k in EXT+CAR},E={k:round(float(XE.loc[i,k]),3) for k in EXT+CAR})
R['peers']=[pr(i) for i in peers]
oy=list(XL.sort_values('Y',ascending=False).index); j=oy.index(TR); R['cargo_peers']=[pr(i) for i in oy[j-6:j+7]]
cp_=oy[j-6:j+7]; R['cargo_peers_med']={k:round(float(XL.loc[[c for c in cp_ if c!=TR],k].median()),3) for k in EXT+CAR+['C','yog']}
R['yog_dist']=dict(n=int((XL.Y>=0.5).sum()),med=round(float(XL[XL.Y>=0.5].yog.median()),3),tr_rank=int(XL[XL.Y>=0.5].yog.rank(ascending=False)[TR]),tr_rank_E=int(XE[XE.Y>=0.5].yog.rank(ascending=False)[TR]),nE=int((XE.Y>=0.5).sum()),medE=round(float(XE[XE.Y>=0.5].yog.median()),3))
# 4) Kaldıraç
q=PER['L']; Wd=dict(insa=Bw[q['flow']].mean(),sahiplik=Ow[q['stock']].mean()*1e3,bayrak=Fw[q['stock']].mean()*1e3,gemiadami=Kw[2021])
R['world_base']={k:float(v) for k,v in Wd.items()}; R['unit_gain']={k:float(0.25*100/v) for k,v in Wd.items()}      # fiziksel birim başına endeks puanı
lev={}; srt=XL.C.sort_values(ascending=False)
for tgt in (20,15,10):
    need=float(srt.iloc[tgt-1]-XL.C[TR]); lev[tgt]=dict(ulke=XL.ad[srt.index[tgt-1]],C=round(float(srt.iloc[tgt-1]),3),fark=round(need,3),gerek={k:dict(pay=round(float(XL.loc[TR,k]+need*4),3),fiziksel=round(need*4/100*Wd[k],0),kat=round(float((XL.loc[TR,k]+need*4)/XL.loc[TR,k]),2)) for k in CAP})
R['lever']=lev
trF=bo[(bo.BeneficialOwnership==TR)&(bo.Year==2026)].set_index('FlagOfRegistration')['Dead weight tons in thousands']; own26=float(trF['0000']); nat26=float(trF['C00001']); R['tr_own26']=own26; R['tr_nat26']=nat26; fw26=float(Fw[2026])
def scen(dx):
    c=XL.C[TR]+sum(dx.values())/4; v=XL.C.copy(); v[TR]=c; return dict(C=round(float(c),3),sira=int(v.rank(ascending=False,method='min')[TR]),dC=round(float(c-XL.C[TR]),3),yog=round(float(c/XL.Y[TR]),3))
half=(own26-nat26)*0.5/fw26*100; crew15=float(XE.loc[TR,'gemiadami']-XL.loc[TR,'gemiadami'])
R['scen']={'bayrak_yari':dict(dpay=round(half,3),**scen({'bayrak':half})),'bayrak_tam':dict(dpay=round(2*half,3),**scen({'bayrak':2*half})),
 'insa_2kat':dict(dpay=round(float(XL.loc[TR,'insa']),3),**scen({'insa':XL.loc[TR,'insa']})),'insa_5kat':dict(dpay=round(float(4*XL.loc[TR,'insa']),3),**scen({'insa':4*XL.loc[TR,'insa']})),
 'gemiadami_2015':dict(dpay=round(crew15,3),**scen({'gemiadami':crew15})),'sahiplik_1_5kat':dict(dpay=round(float(XL.loc[TR,'sahiplik']*.5),3),**scen({'sahiplik':XL.loc[TR,'sahiplik']*.5})),
 'birlikte':scen({'bayrak':half,'insa':XL.loc[TR,'insa'],'gemiadami':crew15})}
# 5) Yıllık seri
ann=pd.DataFrame({'insa':B.loc[TR]/Bw*100,'sahiplik':(O.loc[TR]/Ow*100).shift(-1),'bayrak':(Fl.loc[TR]/Fw*100).reindex(range(2014,2028)).shift(-1),'geri_donusum':S.loc[TR]/Sw*100,'yuk':T_.loc[TR]/Tw*100,'teu':C_.loc[TR]/Cw*100}).loc[2014:2025]
ann['C3']=ann[['insa','sahiplik','bayrak']].mean(axis=1); ann['Y']=ann[['yuk','teu']].mean(axis=1)
ann.round(4).to_csv(O_+'turkiye_yillik.csv'); R['ann']={int(y):{k:(None if pd.isna(v) else round(float(v),3)) for k,v in r.items()} for y,r in ann.iterrows()}
R['tr_crew']={int(y):dict(toplam=int(K_.loc[TR,y]),zabit=int(Ko.loc[TR,y]),pay=round(float(K_.loc[TR,y]/Kw[y]*100),3),sira=int(K_[y].rank(ascending=False,method='min')[TR])) for y in (2015,2021)}; R['world_crew']={int(y):int(Kw[y]) for y in (2015,2021)}
R['tr_phys']=dict(build_gt=float(B.loc[TR,[2023,2024,2025]].mean()),build_gt_E=float(B.loc[TR,[2014,2015,2016]].mean()),own_dwt=own26*1e3,own_dwt_17=float(O.loc[TR,2017])*1e3,flag_dwt=float(Fl.loc[TR,2026])*1e3,flag_dwt_17=float(Fl.loc[TR,2017])*1e3,trade_t=float(T_.loc[TR,[2022,2023,2024]].mean())*1e3,trade_t_E=float(T_.loc[TR,[2014,2015,2016]].mean())*1e3,teu=float(C_.loc[TR,[2022,2023,2024]].mean()),teu_E=float(C_.loc[TR,[2014,2015,2016]].mean()))
C3=['insa','sahiplik','bayrak']
R['no_crew']={p:dict(C3=round(float(X.loc[TR,C3].mean()),3),yog3=round(float(X.loc[TR,C3].mean()/X.Y[TR]),3),sira=int(X[C3].mean(axis=1).rank(ascending=False,method='min')[TR])) for p,X in (('E',XE),('L',XL))}
R['crew_change_big']=[dict(ad=XL.ad[i],E=int(K_.loc[i,2015]),L=int(K_.loc[i,2021])) for i in ['826','208','300','643','724','156','792','608','356']]
R['tr_share_of_C']={k:round(float(XL.loc[TR,k]/4/XL.C[TR]*100),1) for k in CAP}; R['tr_share_of_C5']={k:round(float(XL.loc[TR,k]/5/XL.C5[TR]*100),1) for k in EXT}
json.dump(R,open(O_+'results.json','w'),ensure_ascii=False,indent=1)
assert 'NaN' not in open(O_+'results.json').read()
# ---------------- Şekiller
plt.rcParams.update({'font.size':9,'font.family':'DejaVu Sans','axes.spines.top':False,'axes.spines.right':False}); K1,K2,K3,K4,KG='#1f4e79','#c0392b','#2e8b57','#b8860b','#888888'
LB={'insa':'İnşa','sahiplik':'Sahiplik','bayrak':'Bayrak','gemiadami':'Gemiadamı','geri_donusum':'Geri dönüşüm','yuk':'Deniz yoluyla taşınan yük','teu':'Konteyner elleçleme'}
# Ş1
d=XL[(XL.C>0.02)&(XL.Y>0.02)]; fig,ax=plt.subplots(figsize=(7.2,5.6)); ax.scatter(d.Y,d.C,s=14,c=KG,alpha=.6,lw=0)
lo,hi=0.02,40; ax.plot([lo,hi],[lo,hi],'--',c='k',lw=.7); ax.text(14,19,'C = Y',fontsize=8,rotation=38)
show=['156','410','392','300','430','591','608','702','792','344','584','643','360','208','840','470','826','276','804','704','076','578','784','196','458','036','528','724','380','484','356','050','586','250','756','710','682','818']
from adjustText import adjust_text
ax.set_xscale('log'); ax.set_yscale('log'); ax.set_xlim(lo,hi); ax.set_ylim(lo,hi); tx=[]
for i in show:
    if i in d.index: ax.scatter(d.Y[i],d.C[i],s=42 if i==TR else 18,c=K2 if i==TR else K1,zorder=3); tx.append(ax.text(d.Y[i],d.C[i],d.ad[i],fontsize=8.5 if i==TR else 6.8,fontweight='bold' if i==TR else 'normal',color=K2 if i==TR else 'k'))
adjust_text(tx,ax=ax,arrowprops=dict(arrowstyle='-',color='#999999',lw=.4),expand=(1.15,1.3))
ax.set_xscale('log'); ax.set_yscale('log'); ax.set_xlim(lo,hi); ax.set_ylim(lo,hi); ax.set_xlabel('Yük tabanı endeksi Y (dünya payı, %)'); ax.set_ylabel('Çekirdek denizcilik yetkinliği endeksi C (dünya payı, %)'); fig.tight_layout(); fig.savefig(F_+'Figure_1.png',dpi=300); plt.close()
# Ş2
m=MC.head(30).iloc[::-1]; fig,ax=plt.subplots(figsize=(7.2,7)); yy=np.arange(len(m))
ax.hlines(yy,m.p05,m.p95,color=[K2 if i==TR else K1 for i in m.index],lw=2.2,alpha=.75); ax.scatter(m.med,yy,c='w',edgecolor='k',s=26,zorder=3,label='Ortanca sıra'); ax.scatter(m.r0,yy,marker='|',c='k',s=70,zorder=4,label='Temel endeks sırası'); ax.scatter(m.r5,yy,marker='x',c=K4,s=24,zorder=4,label='Geri dönüşüm dahil sıra')
ax.set_yticks(yy); ax.set_yticklabels(m.ad); [t.set_color(K2) or t.set_fontweight('bold') for t,i in zip(ax.get_yticklabels(),m.index) if i==TR]
ax.set_xlabel('Sıra (1 = en yüksek); çizgi: benzetimlerin %5–%95 aralığı'); ax.set_xlim(0,86); ax.legend(frameon=False,loc='lower right'); fig.tight_layout(); fig.savefig(F_+'Figure_2.png',dpi=300); plt.close()
# Ş3
fig,axs=plt.subplots(1,2,figsize=(7.4,3.4)); a=axs[0]
for k,c in zip(['sahiplik','bayrak','insa'],[K1,K2,K3]): a.plot(ann.index,ann[k],marker='o',ms=3,c=c,label=LB[k])
a.scatter([2015,2021],[R['tr_crew'][2015]['pay'],R['tr_crew'][2021]['pay']],marker='s',c=K4,s=28,label='Gemiadamı (2015, 2021)',zorder=3); a.set_ylabel('Dünya payı (%)'); a.set_title('(a) Yetkinlik göstergeleri',fontsize=9); a.legend(frameon=False,fontsize=7); a.set_ylim(0,2.6)
a=axs[1]
for k,c in zip(['yuk','teu'],[K1,K3]): a.plot(ann.index,ann[k],marker='o',ms=3,c=c,label=LB[k])
a.set_title('(b) Yük tabanı göstergeleri',fontsize=9); a.legend(frameon=False,fontsize=7); a.set_ylim(0,2.6); fig.tight_layout(); fig.savefig(F_+'Figure_3.png',dpi=300); plt.close()
# Ş4
P=pd.DataFrame(R['peers']).set_index('ad'); Kc=pd.DataFrame(list(P.kat3),index=P.index); Kc['top']=Kc.sum(axis=1); Kc=Kc.sort_values('top'); fig,ax=plt.subplots(figsize=(7.2,4.8))
for sign in (1,-1):
    base=np.zeros(len(Kc))
    for k,c in zip(['insa','sahiplik','bayrak'],[K3,K1,K2]):
        v=Kc[k].where(Kc[k]*sign>0,0).values; ax.barh(Kc.index,v,left=base,color=c,label=LB[k] if sign==1 else None,height=.65); base+=v
ax.scatter(Kc.top,np.arange(len(Kc)),c='k',s=14,zorder=3,label='Toplam değişim'); ax.axvline(0,c='k',lw=.6); ax.set_xlabel('Üç göstergeli endeksteki değişim, 2014–2016 → 2023–2025 (yüzde puan)'); ax.legend(frameon=False,fontsize=8,loc='upper left'); [t.set_color(K2) or t.set_fontweight('bold') for t in ax.get_yticklabels() if t.get_text()=='Türkiye']; fig.tight_layout(); fig.savefig(F_+'Figure_4.png',dpi=300); plt.close()
# Ş5
g=R['lever'][15]['gerek']; fig,ax=plt.subplots(figsize=(6.4,3.2)); ks=CAP; vals=[g[k]['kat'] for k in ks]; ax.bar([LB[k] for k in ks],vals,color=[K3,K1,K2,K4],width=.6)
for x,v in enumerate(vals): ax.text(x,v+.15,f'{v:.1f} kat'.replace('.',','),ha='center',fontsize=8)
ax.axhline(1,c='k',lw=.6,ls='--'); ax.set_ylabel('Bugünkü payın katı'); fig.tight_layout(); fig.savefig(F_+'Figure_5.png',dpi=300); plt.close()
t=R['tr']; print({k:t['E'][k] for k in EXT+CAR+['C','C5','Y','yog','r_C','r_C5','r_Y']}); print({k:t['L'][k] for k in EXT+CAR+['C','C5','Y','yog','r_C','r_C5','r_Y']+['r_'+x for x in EXT+CAR]})
print('MC',R['mc_tr'],R['mc_tr_by_opt'],R['mc_tr_S1']); print(R['alt_rank_tr'],R['alt_rank_tr_E']); print('REG',R['reg']); print('LSCI',R['lsci_corr'],R['yog_dist'],R['cargo_peers_med'])
print(pd.DataFrame(R['top30'])[['ad','C','r_C','C_E','rC_E','r_C5','Y','r_Y','yog']+CAP].to_string())
print(pd.DataFrame(R['mc'])[['ad','r0','r5','med','p05','p95']].to_string())
print(pd.DataFrame(R['cargo_peers'])[['ad','C_L','Y_L','yog_L','yog_E','rC_L','rY_L']+CAP].to_string()); print(pd.DataFrame(R['peers'])[['ad','C_E','C_L','rC_E','rC_L','rC5_L','Y_E','Y_L','yog_E','yog_L','art','kat']].to_string())
print(json.dumps(R['lever'],ensure_ascii=False)); print(R['scen']); print(R['res_top']); print(R['res_bot']); print(R['res_big']); print(R['tr_crew'],R['world_crew'],R['hhi'],R['top5share'],R['unit_gain'],R['tr_phys'],R['tr_own26'],R['tr_nat26'],R['ann'])

# ------------------------------------------------------------------------------
# ## Şekiller / Figures
# ------------------------------------------------------------------------------
print('Sekiller / figures:', sorted(glob.glob('fig/*.png')))