def METHOD_NAME(samples, counts, counts_out, prefix, METHOD_NAME, stats, corr):
    # read list of file names
    dat = pd.read_csv(samples, sep='\t', header=None, names=['psample', 'path'])
    # genes-samples matrix
    dfs = [read_counts_tpm_megsap(row.path, row.psample) for _, row in dat.iterrows()]
    d_tpm = dfs[0].join(dfs[1:], how='outer')
    if METHOD_NAME:
        METHOD_NAME.write('#')
        d_tpm.to_csv(METHOD_NAME, sep='\t')
    # statistics per gene
    df_stats = pd.concat([
        d_tpm.mean(axis=1),
        np.log2(d_tpm+1).mean(axis=1),
        np.log2(d_tpm+1).std(axis=1) ],
        axis=1,
        keys=["cohort_mean", "cohort_meanlog2", "cohort_sdlog2"]
        )
    if stats:
        stats.write('#')
        df_stats.to_csv(stats, sep='\t')
    if counts:
        counts_tbl =  pd.read_csv(counts, sep="\t", index_col=0)
        # calculate ratios per gene for sample
        expr = counts_tbl[['tpm']].join(df_stats)
        expr['log2tpm'] = np.log2(expr['tpm']+1)
        expr["log2fc"] = np.log2(expr["tpm"] + 1) - np.log2(expr["cohort_mean"] + 1)
        expr['zscore'] = (np.log2(expr["tpm"] + 1) - expr['cohort_meanlog2'])/expr['cohort_sdlog2']
        expr['pval'] = expr['zscore'].apply(lambda z: 2*norm.cdf(-abs(z)))
        # annotated count file
        if counts_out:
            counts_annot = counts_tbl.join(expr.drop('tpm', axis=1).rename(columns=lambda s: prefix + str(s)),
                how="left")
            counts_annot.to_csv(counts_out, sep='\t', na_rep='n/a')
        if corr:
            value = expr[['log2tpm', 'cohort_meanlog2']].query('log2tpm>0 | cohort_meanlog2>0').corr(method='spearman').iloc[0,1]
            corr.write(f'{value}\n')