plynk 0.1.0

Easy command over the PLINK software directly from Python

README

BETA

This is still undergoing development - Use at your own risk! Notably, there are:

1. NO unittests (will probably be delivered soon). Basic unittesting only. Stuff might break!
2. NO proper testing of all available plink commands.
3. PlinkFileReader.load_file can only read ONE file - Use PlinkFileReader.read_file with a specific file path. args and kwargs follow pandas.read_csv.

Installation

1. Set up a new conda environment (assuing name "my_env").
   1. conda create --name my_env python=3.12 -y
   2. conda activate my_env
2. Install plynk.
   1. This always works: pip install plynk/.
   2. This is faster but only works if you installed uv: uv pip install plynk/. Install uv with pip install uv.

Using PLINK in python to process genetic data

Plink is an abstraction which allows you to use Plink commands easily from Python, with automatic management of path variables so you maintain full control.

Import the Plink class

from library.converters.plink import Plink

Plink is an abstraction which allows you to use the Plink software commands easily from Python, with automatic management of path variables so you maintain full control.

Plink let's you choose three optional parameters:

1. plink_binary_path: This is the path to where you have your plink binary available. If you leave it empty, you must have plink globally installed.
2. plink_prefix: This is the path to the folder where your output data lives, starting with the initial files to process. If you leave it blank, you won't be able to use the file reader. Remember to always create/keep a backup of your source data.
3. encoding: Default is utf-8. If your data is encoded in a different format, consider setting this here.

In the following example, we want to work with a local version of plink which should be at "binaries/plink" in our working directory. Our data is stored in the folder "data_src/demo_data".

# Constants
PLINK_BINARY_PATH = "binaries/plink"
PLINK_PREFIX = "data_src/demo_data/demo_01" # folder: data_scr/demo_data, file_prefix: demo_01

plink = Plink(
    plink_binary_path=PLINK_BINARY_PATH,
    plink_prefix=PLINK_PREFIX,
)

try:
    print(plink.info)
except RuntimeError:
    print("Plink is not yet available. Let's download it!")

Plink is not yet available. Let's download it!

Get the plink binary

All versions of Plink are available here: https://www.cog-genomics.org/plink2/

Plink can download this data for us so that we do not need to take care about it.

plink.download_binaries()

PosixPath('/home/achwalt/Development/plynk/binaries/plink')

plink.info

{'name': 'PLINK',
 'version': 'v1.90b7.2',
 'architecture': '64-bit',
 'release_date': datetime.datetime(2023, 12, 11, 0, 0)}

Working with data

Now that plink is available, we can work with it natively from python. We can use either

1. the exposed Path object or
2. the plink file reader

to see what data we have in our data folder to work with. Since the plink prefix contains already both the folder and the file's prefix, we need to use it's parent!

files_available = [file_name for file_name in plink.plink_prefix.parent.iterdir()]
files_available

[PosixPath('/home/achwalt/Development/plynk/data_src/demo_data/demo_01.ped'),
 PosixPath('/home/achwalt/Development/plynk/data_src/demo_data/demo_01.map')]

plink.file_reader.find_files(plink.plink_prefix)

[PosixPath('/home/achwalt/Development/plynk/data_src/demo_data/demo_01.map'),
 PosixPath('/home/achwalt/Development/plynk/data_src/demo_data/demo_01.ped')]

We can inspect the data by using the path:

ped_data = plink.file_reader.read_file(files_available[0])
ped_data

<style scoped> .dataframe tbody tr th:only-of-type { vertical-align: middle; }

.dataframe tbody tr th {
    vertical-align: top;
}

.dataframe thead th {
    text-align: right;
}

</style>

	TYPE_3c64f9b442261a43a88a89606d363421	COW_3c935a5fb77c7462d6535787c846bd24	0	0.1	0.2	-9	A	G	C	A.1	...	C.7482	C.7483	A.37157	A.37158	G.31809	G.31810	G.31811	G.31812	G.31813	A.37159
0	TYPE_3c64f9b442261a43a88a89606d363421	COW_8ffcd41500450a2afc2a7e18740709d1	0	0	0	-9	A	G	C	A	...	C	C	A	A	G	A	G	G	G	G
1	TYPE_3c64f9b442261a43a88a89606d363421	COW_133f937b04029126ee01146a0c1bb594	0	0	0	-9	A	G	C	A	...	C	C	A	A	G	G	G	G	G	G
2	TYPE_3c64f9b442261a43a88a89606d363421	COW_7d49f51b714f3dbd2c25a04254295b5d	0	0	0	-9	A	G	C	A	...	C	C	A	A	G	A	G	A	G	G
3	TYPE_3c64f9b442261a43a88a89606d363421	COW_15e3a8fe1a0172f527dfb8451492c671	0	0	0	-9	A	G	C	A	...	C	C	A	A	G	A	G	G	G	G
4	TYPE_3c64f9b442261a43a88a89606d363421	COW_08020b767145f321d189fcac4f94cd1c	0	0	0	-9	G	G	C	A	...	A	C	A	A	G	A	G	A	G	G
...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...
60	TYPE_3cda26e0c8aedbf662adb42f923ef3ec	COW_63b051a5b94f870a92026fa87043f7a1	0	0	0	-9	G	G	A	A	...	A	C	A	A	G	A	G	G	G	A
61	TYPE_3cda26e0c8aedbf662adb42f923ef3ec	COW_ff04e2afe8a65c7ab4a6cb89d3be51f8	0	0	0	-9	G	G	A	A	...	C	C	A	A	G	A	G	A	G	A
62	TYPE_3cda26e0c8aedbf662adb42f923ef3ec	COW_4a2448cc04232570fee82a7d68d94bda	0	0	0	-9	G	G	A	A	...	C	C	A	A	G	G	G	G	G	G
63	TYPE_3cda26e0c8aedbf662adb42f923ef3ec	COW_5346afe3589590c3fe017de8c4b8cad5	0	0	0	-9	G	G	A	A	...	A	C	A	A	G	G	G	G	G	A
64	TYPE_3cda26e0c8aedbf662adb42f923ef3ec	COW_e507d53129af542a4f49f4fc14299061	0	0	0	-9	G	G	A	A	...	A	A	A	A	A	A	G	G	G	A

65 rows × 76694 columns

If there is only one file per plink file type, we can also just load it with some pre-formatted column names.

ped_data = plink.file_reader.load_data("ped", numerated_colname="Nucleotide")
ped_data

<style scoped> .dataframe tbody tr th:only-of-type { vertical-align: middle; }

.dataframe tbody tr th {
    vertical-align: top;
}

.dataframe thead th {
    text-align: right;
}

</style>

	Family ID	Individual ID	Paternal ID	Maternal ID	Sex	Phenotype	Nucleotide 1	Nucleotide 2	Nucleotide 3	Nucleotide 4	...	Nucleotide 76679	Nucleotide 76680	Nucleotide 76681	Nucleotide 76682	Nucleotide 76683	Nucleotide 76684	Nucleotide 76685	Nucleotide 76686	Nucleotide 76687	Nucleotide 76688
0	TYPE_3c64f9b442261a43a88a89606d363421	COW_3c935a5fb77c7462d6535787c846bd24	0	0	0	-9	A	G	C	A	...	C	C	A	A	G	G	G	G	G	A
1	TYPE_3c64f9b442261a43a88a89606d363421	COW_8ffcd41500450a2afc2a7e18740709d1	0	0	0	-9	A	G	C	A	...	C	C	A	A	G	A	G	G	G	G
2	TYPE_3c64f9b442261a43a88a89606d363421	COW_133f937b04029126ee01146a0c1bb594	0	0	0	-9	A	G	C	A	...	C	C	A	A	G	G	G	G	G	G
3	TYPE_3c64f9b442261a43a88a89606d363421	COW_7d49f51b714f3dbd2c25a04254295b5d	0	0	0	-9	A	G	C	A	...	C	C	A	A	G	A	G	A	G	G
4	TYPE_3c64f9b442261a43a88a89606d363421	COW_15e3a8fe1a0172f527dfb8451492c671	0	0	0	-9	A	G	C	A	...	C	C	A	A	G	A	G	G	G	G
...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...
61	TYPE_3cda26e0c8aedbf662adb42f923ef3ec	COW_63b051a5b94f870a92026fa87043f7a1	0	0	0	-9	G	G	A	A	...	A	C	A	A	G	A	G	G	G	A
62	TYPE_3cda26e0c8aedbf662adb42f923ef3ec	COW_ff04e2afe8a65c7ab4a6cb89d3be51f8	0	0	0	-9	G	G	A	A	...	C	C	A	A	G	A	G	A	G	A
63	TYPE_3cda26e0c8aedbf662adb42f923ef3ec	COW_4a2448cc04232570fee82a7d68d94bda	0	0	0	-9	G	G	A	A	...	C	C	A	A	G	G	G	G	G	G
64	TYPE_3cda26e0c8aedbf662adb42f923ef3ec	COW_5346afe3589590c3fe017de8c4b8cad5	0	0	0	-9	G	G	A	A	...	A	C	A	A	G	G	G	G	G	A
65	TYPE_3cda26e0c8aedbf662adb42f923ef3ec	COW_e507d53129af542a4f49f4fc14299061	0	0	0	-9	G	G	A	A	...	A	A	A	A	A	A	G	G	G	A

66 rows × 76694 columns

map_data = plink.file_reader.load_data("map")
map_data

<style scoped> .dataframe tbody tr th:only-of-type { vertical-align: middle; }

.dataframe tbody tr th {
    vertical-align: top;
}

.dataframe thead th {
    text-align: right;
}

</style>

	Chromosome	Marker ID	Genetic Distance	Position Base Pairs
0	1	COW_bee91d0199f55f7ceac640dd11d140c0	0	135098
1	1	COW_1913080647293a3848ad937ca9653aec	0	149772
2	1	COW_2eb599df7f168d419f8c2757da552295	0	163995
3	1	COW_96317923c4b14b72dc44ebbdfa25d279	0	183040
4	1	COW_81d964d5d890f1b262d1cb2d8e784523	0	267940
...	...	...	...	...
38339	29	COW_7319e9361f794ba225c89a83edfd620b	0	51038808
38340	29	COW_e0c1b1b57e7e6d7b6c83573badbaa991	0	51042316
38341	29	COW_a7f1212ead02b209c58c1766a914b01a	0	51152056
38342	29	COW_5b06fd46612881ff9e7d2604de917a46	0	51358815
38343	29	COW_34d9e4d41e8d8c8b8cadc8505a4f09d6	0	51484561

38344 rows × 4 columns

Processing the data for further analysis

Plink enables you to process your data in order to make it easy to analyse. While Plink itself only exposes a command line API, you must still maintain all used instructions for your analysis in order to be scientific: Research must be reproducible.

This is why this package was created: It allows you to do your analysis in Python while also commanding Plink natively from within Python.

Plink allows you to use pythonic arguments and keyword arguments, which will then be translated into command line instructions.

# Encouraged: Use a pythonic way to declare args and kwargs
output = plink.run(
    make_bed=True,
    file=plink.plink_prefix,
    out=plink.plink_prefix,
    chr_set=34,
    chr="1-29",
    geno=0.1,
    mind=0.1,
    maf=0.05
)

# Supports suppling arguments as strings directly, such as '--make-bed'
output = plink.run(
    "--make-bed",
    file=plink.plink_prefix,
    out=plink.plink_prefix,
    chr_set=34,
    chr="1-29",
    geno=0.1,
    mind=0.1,
    maf=0.05
)

Type Safety for Plink keywords

Plink will prefer type-safe conversion where possible. But you can indeed provide basic arguments such as --make-bed directly as well.

Currently as subsidiary method converts plink commands unsafely to not block you from working due to a missing feature. That's why it is discouraged to provide arguments as strings.

cmd_used = plink.make_cmd(
    "--make-bed",
    file=plink.plink_prefix,
    out=plink.plink_prefix,
    chr_set=34,
    chr="1-29",
    geno=0.1,
    mind=0.1,
    maf=0.05
)
cmd_used

[PosixPath('/home/achwalt/Development/plynk/binaries/plink'),
 '--make-bed',
 <PlinkKeyword.FILE: '--file'>,
 '/home/achwalt/Development/plynk/data_src/demo_data/demo_01',
 <PlinkKeyword.OUT: '--out'>,
 '/home/achwalt/Development/plynk/data_src/demo_data/demo_01',
 <PlinkKeyword.CHR_SET: '--chr-set'>,
 '34',
 <PlinkKeyword.CHR: '--chr'>,
 '1-29',
 <PlinkKeyword.GENO: '--geno'>,
 '0.1',
 <PlinkKeyword.MIND: '--mind'>,
 '0.1',
 <PlinkKeyword.MAF: '--maf'>,
 '0.05']

[str(cmd) for cmd in cmd_used]

['/home/achwalt/Development/plynk/binaries/plink',
 '--make-bed',
 '--file',
 '/home/achwalt/Development/plynk/data_src/demo_data/demo_01',
 '--out',
 '/home/achwalt/Development/plynk/data_src/demo_data/demo_01',
 '--chr-set',
 '34',
 '--chr',
 '1-29',
 '--geno',
 '0.1',
 '--mind',
 '0.1',
 '--maf',
 '0.05']

from pathlib import Path

plink.make_cmd(
    "--make-bed",
    "--file",
    'data_src/demo_data/demo_01', # ❌ No Path safety will be applied!
    "--out",
    Path('data_src/demo_data/demo_01'), # ✅ Path safety will be applied.
    chr_set=34,
    chr="1-29",
    geno=0.1,
    mind=0.1,
    maf=0.05
)

[PosixPath('/home/achwalt/Development/plynk/binaries/plink'),
 '--make-bed',
 '--file',
 'data_src/demo_data/demo_01',
 '--out',
 PosixPath('/home/achwalt/Development/plynk/data_src/demo_data/demo_01'),
 <PlinkKeyword.CHR_SET: '--chr-set'>,
 '34',
 <PlinkKeyword.CHR: '--chr'>,
 '1-29',
 <PlinkKeyword.GENO: '--geno'>,
 '0.1',
 <PlinkKeyword.MIND: '--mind'>,
 '0.1',
 <PlinkKeyword.MAF: '--maf'>,
 '0.05']

You can also run a command created by make_cmd

cmd = plink.make_cmd(
    make_bed=True,
    file=plink.plink_prefix,
    out=plink.plink_prefix,
    chr_set=34,
    chr="1-29",
    geno=0.1,
    mind=0.1,
    maf=0.05
)
output = plink.run_cmd(cmd)

# Let's see what files we have now
plink.file_reader.find_files(plink.plink_prefix)

[PosixPath('/home/achwalt/Development/plynk/data_src/demo_data/demo_01.bed'),
 PosixPath('/home/achwalt/Development/plynk/data_src/demo_data/demo_01.bim'),
 PosixPath('/home/achwalt/Development/plynk/data_src/demo_data/demo_01.fam'),
 PosixPath('/home/achwalt/Development/plynk/data_src/demo_data/demo_01.log'),
 PosixPath('/home/achwalt/Development/plynk/data_src/demo_data/demo_01.map'),
 PosixPath('/home/achwalt/Development/plynk/data_src/demo_data/demo_01.nosex'),
 PosixPath('/home/achwalt/Development/plynk/data_src/demo_data/demo_01.ped')]

## More Examples
cmd = plink.make_cmd(file=plink.plink_prefix, out=plink.plink_prefix, cow=True, het=True)
output = plink.run_cmd(cmd)

plink.file_reader.load_data("het")

<style scoped> .dataframe tbody tr th:only-of-type { vertical-align: middle; }

.dataframe tbody tr th {
    vertical-align: top;
}

.dataframe thead th {
    text-align: right;
}

</style>

	FID	IID	O(HOM)	E(HOM)	N(NM)	F
0	TYPE_3c64f9b442261a43a88a89606d363421	COW_3c935a5fb77c7462d6535787c846bd24	22120	22280.0	38015	-0.010110
1	TYPE_3c64f9b442261a43a88a89606d363421	COW_8ffcd41500450a2afc2a7e18740709d1	22102	22280.0	38017	-0.011320
2	TYPE_3c64f9b442261a43a88a89606d363421	COW_133f937b04029126ee01146a0c1bb594	22005	22280.0	38017	-0.017490
3	TYPE_3c64f9b442261a43a88a89606d363421	COW_7d49f51b714f3dbd2c25a04254295b5d	22194	22280.0	38018	-0.005507
4	TYPE_3c64f9b442261a43a88a89606d363421	COW_15e3a8fe1a0172f527dfb8451492c671	21689	22280.0	38018	-0.037600
...	...	...	...	...	...	...
61	TYPE_3cda26e0c8aedbf662adb42f923ef3ec	COW_63b051a5b94f870a92026fa87043f7a1	23399	22280.0	38018	0.071060
62	TYPE_3cda26e0c8aedbf662adb42f923ef3ec	COW_ff04e2afe8a65c7ab4a6cb89d3be51f8	23515	22280.0	38017	0.078470
63	TYPE_3cda26e0c8aedbf662adb42f923ef3ec	COW_4a2448cc04232570fee82a7d68d94bda	23970	22280.0	38018	0.107300
64	TYPE_3cda26e0c8aedbf662adb42f923ef3ec	COW_5346afe3589590c3fe017de8c4b8cad5	23874	22280.0	38018	0.101200
65	TYPE_3cda26e0c8aedbf662adb42f923ef3ec	COW_e507d53129af542a4f49f4fc14299061	25778	22280.0	38018	0.222200

66 rows × 6 columns