sempwn / LF-model

Individual-based stochastic simulation of Lymphatic Filariasis infection. The accompanying paper for this model can be found here: http://www.parasitesandvectors.com/content/8/1/547
MIT License
0 stars 0 forks source link

Working example? #1

Open sdwfrost opened 2 years ago

sdwfrost commented 2 years ago

I've tried to run the code after compiling as follows - any help would be appreciated.

Windows:

g++ -O3 class1.cpp modelclasses.cpp -I . -o lf -lgsl -lopenblas

Ubuntu:

g++ -O3 class1.cpp modelclasses.cpp -I . -o lf -lgsl -lgslblas

Running without command line options gives the following output and error:

reading line :1.0 : riskMu1
token :1.0
reading line :1.0 : riskMu2
token :1.0
reading line :1.0 : riskMu3
token :1.0
reading line :0.065 : shapeRisk //shape parameter for bite-risk distribution (0.1/0.065)
token :0.065
reading line :0.0104 : mu //death rate of worms
token :0.0104
reading line :0.0 : theta //0.001 //immune system response parameter. 0.112
token :0.0
reading line :0.1 :  gamma //mf death rate
token :0.1
reading line :0.58 : alpha //mf birth rate per fertile worm per 20 uL of blood.
token :0.58
reading line :10 : lbda //number of bites per mosquito per month.
token :10
reading line :9.0 : v_to_h  //vector to host ratio (39.,60.,120.)
token :9.0
reading line :4.395 :  kappas1 //vector uptake and development anophelene
token :4.395
reading line :0.055 : r1 //vector uptake and development anophelene
token :0.055
reading line :0.00167 : tau  //death rate of population
token :0.00167
reading line :0.0 : z  //waning immunity
token :0.0
reading line :0.0 : nu  //poly-monogamy parameter
token :0.0
reading line :0.0 : L3  //larvae density.
token :0.0
reading line :0.37 : g //Proportion of mosquitoes which pick up infection when biting an infected host
token :0.37
reading line :5.0 : sig //death rate of mosquitos
token :5.0
reading line :0.414 : psi1 //Proportion of L3 leaving mosquito per bite
token :0.414
reading line :0.32 : psi2 //Proportion of L3 leaving mosquito that enter host
token :0.32
reading line :1.0 : dt //time spacing (months)
token :1.0
reading line :1.0 : lbdaR //use of bed-net leading to reduction in bite rate
token :1.0
reading line :1.0 : v_to_hR //use of residual-spraying leading to reduction in v_to_h
token :1.0
reading line :7 : nMDA //number of rounds of MDA
token :7
reading line :12 : mdaFreq //frequency of MDA (months)
token :12
reading line :0.65 : covMDA //coverage of MDA
token :0.65
reading line :0.00275 : s2 //probability of L3 developing into adult worm.
token :0.00275
reading line :0.95 : chi //proportion of mf removed for a single MDA round.
token :0.95
reading line :0.55 : tau //proportion of worms permanently sterilised for a single MDA round. (0.55)
token :0.55
reading line :0.999 : rho //proportion of systematic non-compliance 0- none 1- all.
token :0.999
reading line :0 : mosquitoSpecies // 0 - Anopheles facilitation squared, 1 - Culex limitation linear.
token :0
reading line :0.0 : rhoBU //correlation between bite risk and systematic non-compliance.
token :0.0
reading line :0 : aWol //using doxycycline in intervention 0- not used, 1- is used.
token :0
reading line :5.0 : sigR //new mortality rate of mosquitoes during vector intervention.
token :5.0
reading line :0.0 : covN //coverage of bed nets.
token :0.0
reading line :0.99 : sysCompN //systematic non-compliance of bed nets. set to near one.
token :0.99
reading line :0.0 : rhoCN //correlation between receiving chemotherapy and use of bed nets.
token :0.0
reading line :0 : IDAControl //if 1 then programme switches to IDA after five rounds of standard MDA defined with chi and tau.
token :0
random u , b and n-0.733634 0.764271 inf
random gamma example: 0.0671858
random poisson example: 102
random gaussian example: 0.89849
Initialising array...
gsl: gammainv.c:111: ERROR: inverse failed to converge
Default GSL error handler invoked.
Aborted
sempwn commented 2 years ago

It's been a while since I've looked over this code. Try compiling with the following MakeFile and calling using the below Python method:

OBJS := timer.o randomseed.o

# Provide any options and flags for the compiler here
CFLAGS    = -Wall -O3

# Set name of final executable file here
EXENAME = timer

# If you are linking to external libraries, provide them here
LIBS  =  -lgsl -lgslcblas -lfftw3 -lm

# If there are libraries with non-standard paths, provide the paths here
LIBDIR = -L$(HOME)/lib

# If there are header files with non-standard paths, provide the paths here
INCLDIR = -I. -I$(HOME)/include

default:$(EXENAME)

# link to generate executable file
$(EXENAME): $(OBJS)
    g++ $(LIBDIR) $(OBJS) $(LIBS) -o $(EXENAME)

# pull in dependency info for *existing* .o files
-include $(OBJS:.o=.d)

# compile and generate dependency info
%.o: %.cpp
    gcc -c $(CFLAGS) $(INCLDIR) $*.cpp -o $*.o
    gcc -MM $(CFLAGS) $(INCLDIR) $*.cpp > $*.d

# remove compilation products
clean:
    rm -f $(EXENAME) *.o *.d
def runSimNoInput(filename='test3.txt',finalDistribution=True,mfDistribution=False,allDistribution=False,index=0,model_path=''):
    ''' use for parrellisation as does not change the parameter file. '''
    '''run model '''
    if finalDistribution:
        args = str(yrs) + " " + str(n) + " " + str(filename) + " " + str(index)
    else:
        args = str(120.0) + " " + str(n) + " " + str(filename) + " " + str(index)
    cmd = model_path+'timer'
    ccmd = cmd + " " + args
    subprocess.call(ccmd,shell=True)
    '''import data '''
    outtm = np.loadtxt(model_path + filename)

    wp = (outtm[0::6,:]+outtm[0::6,:])>0
    w = outtm[0::6,:]+outtm[0::6,:]
    mfp = outtm[2::6,:]>1.0
    mf = outtm[2::6,:]
    ages = outtm[-2,:]/12.0
    if finalDistribution:
        wp = wp[-1,:]
        mfp = mfp[-1,:]
        mf = mf[-1,:]
    os.remove(filename)
    ts = np.concatenate((np.linspace(0,100,100),np.linspace(100,120,120)))
    if mfDistribution:
        return ts, mfp, wp, ages, mf
    elif allDistribution:
        return ts, mfp, wp, ages, mf, w
    else:
        return ts, mfp, wp, ages
sdwfrost commented 2 years ago

Neither the Makefile or the Python function run out of the box, but certainly calling the program with 4 arguments helps, but the error is stochastic.

Compile:

g++ -O3 class1.cpp modelclasses.cpp -I . -o lf -lgsl -lgslcblas -lm

Run:

./lf 120 1 test.out 1

Working output:

four arguments imported. 
reading line :1.0 : riskMu1 
token :1.0 
reading line :1.0 : riskMu2 
token :1.0 
reading line :1.0 : riskMu3 
token :1.0 
reading line :0.065 : shapeRisk //shape parameter for bite-risk distribution (0.1/0.065)
token :0.065 
reading line :0.0104 : mu //death rate of worms
token :0.0104 
reading line :0.0 : theta //0.001 //immune system response parameter. 0.112
token :0.0 
reading line :0.1 :  gamma //mf death rate
token :0.1 
reading line :0.58 : alpha //mf birth rate per fertile worm per 20 uL of blood.
token :0.58 
reading line :10 : lbda //number of bites per mosquito per month.
token :10 
reading line :9.0 : v_to_h  //vector to host ratio (39.,60.,120.) 
token :9.0 
reading line :4.395 :  kappas1 //vector uptake and development anophelene
token :4.395 
reading line :0.055 : r1 //vector uptake and development anophelene
token :0.055 
reading line :0.00167 : tau  //death rate of population
token :0.00167 
reading line :0.0 : z  //waning immunity
token :0.0 
reading line :0.0 : nu  //poly-monogamy parameter       
token :0.0 
reading line :0.0 : L3  //larvae density.
token :0.0 
reading line :0.37 : g //Proportion of mosquitoes which pick up infection when biting an infected host
token :0.37 
reading line :5.0 : sig //death rate of mosquitos
token :5.0 
reading line :0.414 : psi1 //Proportion of L3 leaving mosquito per bite
token :0.414 
reading line :0.32 : psi2 //Proportion of L3 leaving mosquito that enter host
token :0.32 
reading line :1.0 : dt //time spacing (months) 
token :1.0 
reading line :1.0 : lbdaR //use of bed-net leading to reduction in bite rate
token :1.0 
reading line :1.0 : v_to_hR //use of residual-spraying leading to reduction in v_to_h 
token :1.0 
reading line :7 : nMDA //number of rounds of MDA
token :7 
reading line :12 : mdaFreq //frequency of MDA (months)
token :12 
reading line :0.65 : covMDA //coverage of MDA
token :0.65 
reading line :0.00275 : s2 //probability of L3 developing into adult worm.
token :0.00275 
reading line :0.95 : chi //proportion of mf removed for a single MDA round.
token :0.95 
reading line :0.55 : tau //proportion of worms permanently sterilised for a single MDA round. (0.55)
token :0.55 
reading line :0.999 : rho //proportion of systematic non-compliance 0- none 1- all.
token :0.999 
reading line :0 : mosquitoSpecies // 0 - Anopheles facilitation squared, 1 - Culex limitation linear.
token :0 
reading line :0.0 : rhoBU //correlation between bite risk and systematic non-compliance.
token :0.0 
reading line :0 : aWol //using doxycycline in intervention 0- not used, 1- is used.
token :0 
reading line :5.0 : sigR //new mortality rate of mosquitoes during vector intervention.
token :5.0 
reading line :0.0 : covN //coverage of bed nets.
token :0.0 
reading line :0.99 : sysCompN //systematic non-compliance of bed nets. set to near one.
token :0.99 
reading line :0.0 : rhoCN //correlation between receiving chemotherapy and use of bed nets. 
token :0.0 
reading line :0 : IDAControl //if 1 then programme switches to IDA after five rounds of standard MDA defined with chi and tau.
token :0 
random u , b and n-0.963077 0.730153 inf
random gamma example: 0.486468
random poisson example: 108
random gaussian example: 0.906707
random gamma example inverse method2.39716
Initialising array... 
mosquito species: 0
0----------100
---------bednet event at 1200--
Finished simulation in 0.004536s
write to file... 
destructing random number generator

Broken (different seed):

four arguments imported. 
reading line :1.0 : riskMu1 
token :1.0 
reading line :1.0 : riskMu2 
token :1.0 
reading line :1.0 : riskMu3 
token :1.0 
reading line :0.065 : shapeRisk //shape parameter for bite-risk distribution (0.1/0.065)
token :0.065 
reading line :0.0104 : mu //death rate of worms
token :0.0104 
reading line :0.0 : theta //0.001 //immune system response parameter. 0.112
token :0.0 
reading line :0.1 :  gamma //mf death rate
token :0.1 
reading line :0.58 : alpha //mf birth rate per fertile worm per 20 uL of blood.
token :0.58 
reading line :10 : lbda //number of bites per mosquito per month.
token :10 
reading line :9.0 : v_to_h  //vector to host ratio (39.,60.,120.) 
token :9.0 
reading line :4.395 :  kappas1 //vector uptake and development anophelene
token :4.395 
reading line :0.055 : r1 //vector uptake and development anophelene
token :0.055 
reading line :0.00167 : tau  //death rate of population
token :0.00167 
reading line :0.0 : z  //waning immunity
token :0.0 
reading line :0.0 : nu  //poly-monogamy parameter       
token :0.0 
reading line :0.0 : L3  //larvae density.
token :0.0 
reading line :0.37 : g //Proportion of mosquitoes which pick up infection when biting an infected host
token :0.37 
reading line :5.0 : sig //death rate of mosquitos
token :5.0 
reading line :0.414 : psi1 //Proportion of L3 leaving mosquito per bite
token :0.414 
reading line :0.32 : psi2 //Proportion of L3 leaving mosquito that enter host
token :0.32 
reading line :1.0 : dt //time spacing (months) 
token :1.0 
reading line :1.0 : lbdaR //use of bed-net leading to reduction in bite rate
token :1.0 
reading line :1.0 : v_to_hR //use of residual-spraying leading to reduction in v_to_h 
token :1.0 
reading line :7 : nMDA //number of rounds of MDA
token :7 
reading line :12 : mdaFreq //frequency of MDA (months)
token :12 
reading line :0.65 : covMDA //coverage of MDA
token :0.65 
reading line :0.00275 : s2 //probability of L3 developing into adult worm.
token :0.00275 
reading line :0.95 : chi //proportion of mf removed for a single MDA round.
token :0.95 
reading line :0.55 : tau //proportion of worms permanently sterilised for a single MDA round. (0.55)
token :0.55 
reading line :0.999 : rho //proportion of systematic non-compliance 0- none 1- all.
token :0.999 
reading line :0 : mosquitoSpecies // 0 - Anopheles facilitation squared, 1 - Culex limitation linear.
token :0 
reading line :0.0 : rhoBU //correlation between bite risk and systematic non-compliance.
token :0.0 
reading line :0 : aWol //using doxycycline in intervention 0- not used, 1- is used.
token :0 
reading line :5.0 : sigR //new mortality rate of mosquitoes during vector intervention.
token :5.0 
reading line :0.0 : covN //coverage of bed nets.
token :0.0 
reading line :0.99 : sysCompN //systematic non-compliance of bed nets. set to near one.
token :0.99 
reading line :0.0 : rhoCN //correlation between receiving chemotherapy and use of bed nets. 
token :0.0 
reading line :0 : IDAControl //if 1 then programme switches to IDA after five rounds of standard MDA defined with chi and tau.
token :0 
random u , b and n-1.55928 0.205291 inf
random gamma example: 0.155269
random poisson example: 105
random gaussian example: 1.13217
random gamma example inverse method1.32284
Initialising array... 
gsl: gammainv.c:111: ERROR: inverse failed to converge
Default GSL error handler invoked.
Aborted (core dumped)