How to obtain HomoMating within a set of virtual subpopulations

[ashander]

As the code below demonstrates, running HomoMating with subPops=ALL_AVAIL vs explicitly specifying all virtual subpopulations gives different results. (This is with a single physical subpopulation.)

Is this by design? The docs on argument subPopSize state

If a list of (virtual) subpopulation is specified, the mating scheme will be applied to specific (virtual) subpopulations.

which isn't entirely clear.

If it is by design, is there a way to obtain HomoMating among all individuals in a virtual subpopulation (within the same physical subpopulation)?

Thanks!

from simuPOP import *

def sc(sex):
    if sex == MALE:
        return 'M'
    else:
        return 'F'

for sps in (
    ALL_AVAIL,
    [(ALL_AVAIL, 0), (ALL_AVAIL, 1)],
    [(0, 0), (0, 1)],
    ((ALL_AVAIL, 0), (ALL_AVAIL, 1)),
):
    print('\n\nTrying argument subPops={}'.format(sps))
    try:
        pop = Population([4] * 1, infoFields = ['mark', 'ind_id', 'father_id', 'mother_id'])
        initInfo(pop, values=[0, 1], infoFields= ['mark'])
        initSex(pop, sex=[MALE, FEMALE, MALE, FEMALE])
        pop.setVirtualSplitter(InfoSplitter('mark', values=[0, 1]))
        pop.numVirtualSubPop()
        print(' Ind: {} ({}) {}'.format(' ', 'sex', 'vsp'))
        for i in pop.allIndividuals(subPops=sps):
            print(' Ind: {} ({}) {}'.format(i.ind_id, sc(i.sex()), i.mark))
        pop.evolve(initOps = [IdTagger()],
                   matingScheme = HeteroMating(
                       [
                           HomoMating(
                               subPops=sps,
                               chooser = RandomParentsChooser(),
                               generator = OffspringGenerator(
                                   ops = [
                                       IdTagger(),
                                       PedigreeTagger(),
                                       MendelianGenoTransmitter(),
                                   ],
                                   numOffspring=2,
                               )
                           )
                       ]
                   ),
                   gen=1)
    except Exception as e:
        print(str(e))

Output

Trying argument subPops=True
 Ind:   (sex) vsp
 Ind: 0.0 (M) 0.0
 Ind: 0.0 (F) 1.0
 Ind: 0.0 (M) 0.0
 Ind: 0.0 (F) 1.0

Trying argument subPops=[(True, 0), (True, 1)]
 Ind:   (sex) vsp
 Ind: 0.0 (M) 0.0
 Ind: 0.0 (M) 0.0
 Ind: 0.0 (F) 1.0
 Ind: 0.0 (F) 1.0
Mating scheme 0 in subpopulation 0 failed to produce 2 offspring.
RandomParentsChooser fails because there is no female individual in a subpopulation 

Trying argument subPops=[(0, 0), (0, 1)]
 Ind:   (sex) vsp
 Ind: 0.0 (M) 0.0
 Ind: 0.0 (M) 0.0
 Ind: 0.0 (F) 1.0
 Ind: 0.0 (F) 1.0
Mating scheme 0 in subpopulation 0 failed to produce 2 offspring.
RandomParentsChooser fails because there is no female individual in a subpopulation 

Trying argument subPops=((True, 0), (True, 1))
 Ind:   (sex) vsp
 Ind: 0.0 (M) 0.0
 Ind: 0.0 (M) 0.0
 Ind: 0.0 (F) 1.0
 Ind: 0.0 (F) 1.0
Mating scheme 0 in subpopulation 0 failed to produce 2 offspring.
RandomParentsChooser fails because there is no female individual in a subpopulation 

[petrelharp]

@ashander - you expect (and want) HomoMating to be applied to all the individuals, grouped together, but it seems that HomoMating is being applied separately to each virtual subpopulation?

This is so that eventually we can have mating among everyone in a set of virtual subpops, within each physical subpop.

[BoPeng]

Without checking your code carefully,

• ALL_AVAIL expands to all available subpopulations (not virtual subpopulation), so in the single subpopulation case this is simply [0].
• [(ALL_AVAIL, ALL_AVAIL)] expands to all virtual subpopulations in all subpopulations, so in your case is equivalent to [(0, ALL_AVAIL)], and to [(0, 0), (0,1)].
• [ (ALL_AVAIL, 0), (ALL_AVAIL, 1)] in this case is basically [(0, 0), (0,1)]

So in the first case ALL_AVAIL all individuals in the first subpopulation are involved in mating. In the second case mating happens in two vsps separately. If you happen to define vsp by sex, no asexual mating can happen in these vsps.

[ashander]

Great, that's what I thought was happening based on the output.

But @petrelharp is correct that my ultimate goal is to get mating among everyone set of vsps (call these "mating" vsps) which is a subset of those that exist (call the others "non-mating") instead of separately. Is it true that there not a way to do this via the subPops argument?

If so, would the best way to do it be by specifying fitness =0 in the "non-mating" virtual subpops, or is there some other means?

Thanks,

[BoPeng]

If you have two vsps and would like to mate in the second, use subPops=[(0,1)].

[petrelharp]

Say there are three virtual subpopulations; the first one does not mate, but the second two can mate with each other. So if VSP 0 : 3 males, 3 females VSP 1: 2 males, 4 females VSP 2: 5 males then what @ashander intended would be to call subPops=[(0,1),(0,2)] and having mating occur between all 7 males and 4 females, together.

However: the point of VSPs is they define who mates together - @ashander, do we actually need to do this?

[BoPeng]

Mating happens within subpopulation, virtual or not, so subPops=[(0,1),(0,2)] will not cut it if you would like to mate between all individuals in these two vsps. In this case, you can either redefine your vsps or define a new vsp to combine (0,1) and (0,2) (remember vsps are "definitions" and can overlap with each other). A CombinedSplitter can be really useful here.

[ashander]

An set of VSPs that defines non-mating / mating would do the trick. Thanks @BoPeng

For @petrelharp example: we can define: VSP 3: 7 males, 4 females (contains both (0,1) and (0,2)) and have mating occur only in VSP 3

[ashander]

For completeness, defining a new pair of VSPs for immature and mature individuals works:

pop = Population([4] * 1, infoFields = ['mark', 'ind_id', 'father_id', 'mother_id'])
initInfo(pop, values=[0, 1, 2, 3], infoFields= ['mark'])
initSex(pop, sex=[MALE, FEMALE, MALE, FEMALE])

# some underlying structure that we don't want to affect mating
pop.setVirtualSplitter(InfoSplitter('mark', values=[0, 1, 2, 3]))
num_initial_virtual_subpops = pop.numVirtualSubPop()
# split the population into two VSPs
maturity_splitter = InfoSplitter('mark', cutoff=2)
immature_idx = num_initial_virtual_subpops
mature_idx = num_initial_virtual_subpops + 1
pop.setVirtualSplitter(
    CombinedSplitter([pop.virtualSplitter(),
                      maturity_splitter]))
sps = [(ALL_AVAIL, mature_idx), ]
pop.evolve(initOps = [IdTagger()],
           matingScheme = HeteroMating(
               [
                   HomoMating(
                       subPops=sps,
                       chooser = RandomParentsChooser(),
                       generator = OffspringGenerator(
                           ops = [
                               IdTagger(),
                               PedigreeTagger(),
                               MendelianGenoTransmitter(),
                           ],
                           numOffspring=2,
                       )
                   )
               ]
           ),
           gen=1)